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ONE-DIMENSIONAL  TWO-PHASE  FLOW  IN  COMBUSTION  CHAMBER 
OF  SOLID  PROPELLANT  ROCKET  MOTORS 

CHANG  XIANQI 

ABSTRACT 

'■In  this  paper,  a  numerical  solution  of  basic  equation  for 
one  dimensional  two-phase  nonequilibrium  flow  in  a  combustion 
chamber  of  solid  propellant  rocket  motors  is  discussed  in  detail, 

combust  ion 

the  effect  of  particle  size  on  flow  field  in  chamber  and  pressure 

time  curves  is  analyzed,  and  some  useful  conclusions  are  obtained 

in  comparison  with  results  of  one  dimensional  two-phase  constant 
combustion 

lag  flow  in. chamber.  It  is  useful  for  predicting  pressure-time 
/» 

curves  accurately  and  providing  accurate  boundary  conditions 
for  the  calculation  of  two-phase  flow  through  the  nozzle.  (  H',  •. 


Nomenclature 


/II 


A  -  duct  cross  section  area 
-  Burning  area  of  charge 
Aq  -  Thermal  equivalent  of  work 
Ai  -  Area  of  the  nozzle  throat 
M  -  Mean  molecular  weight  of  gas 
n  -  Pressure  exponent 
P  -  Pressure 
Ps  -  Total  pressure 
b  *  Burning  speed  coefficient 
C*  -  Characteristic  speed  of  propellant 

Ct  -  Particulate  specific  heat 

Cpg  -  Specific  heat  of  gaseous  phase  at  constant  pressure 

g  -  Gravitational  acceleration 
h  -  Enthalpy  per  unit  mass 

Hs  -  Total  enthalpy  in  1  kg  of  two-phase  mixture 
Hsg  -  Total  enthalpy  in  1  kg  of  gas 
Hsp  -  Total  enthalpy  in  1  kg  of  liquid 
k  -  Ratio  of  specific  heat  of  gas 
K  -  Particle  velocity  lag  coefficient 

A, 

K\  -  Surface  -  throat  ratio  J\.  =  - 

At 

1  -  Charge  length 

L-  Particle  temperature  lag  coefficient 
m  -  Mass  flow  rate  of  two-phase  mixture 
rrig  -  Mass  flow  rate  of  gas 
q  -  Heat  flux  of  particles  per  unit  mass 

r  -  Burning  rate 
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rp  -  Radius  of  particles 

Rg  -  Gas  constant  of  gaseous  phase 

S  -  Duct  circumference 
t  -  Time 

T  -  Temperature 

To  ~  Total  temperature  of  charge  head  gaseous  phase 
Ts  -  Total  temperature  of  gaseous  phase 
v -  Velocity 
x  -  Axial  coordinate 
X  -  Particle  resistance  per  unit  mass 
yc  -  Outer  radius  of  charge 
p  -  Density 

pt  -  Density  of  propellant 
771^  -  Mass  flow  rate  of  condensed  phase 
Mg  -  Mach  number  of  gaseous  phase 
Pmp  ‘  Density  of  AL20 3  material 

c  -  Fraction  of  particle  mass  flow  rate,  e  = 

771 

Ay  -  Coefficient  of  thermal  conductivity  of  gaseous  phase 
p.g  -  Coefficient  of  dynamic  viscosity  of  gaseous  phase 

SUBSCRIPTS 

-  Gas  phase 

o  -  Cross  section  of  charge  head 
t-  Cross  section  of  nozzle  throat 
P  -  Condensed  piiase 

L  -  Cross  section  at  the  charge  tip 

i  -  initial  value 
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1.  PREFACE 


Modern  composite  solid  propellants  are  mixed  with  a  certain 
amount  of  aluminum  powder  for  increasing  energy  and  decreasing 
unsteady  burning.  When  the  aluminum  mixed  propellant  burns, 

Al203  particles  are  formed  in  liquid  phase.  The  weight  fraction 
can  be  up  to  30-40%.  Therefore,  in  combustion  chambers  and  jet 
nozzles,  the  products  of  combustion  actually  are  a  mixture  of 
gas  and  liquid. 

There  have  been  a  lot  of  publications  regarding  two-phase 
flow  in  jet  nozzles  of  a  rocket  motor.  Most  of  them  are  devoted 
to  one  dimensional  two-phase  flow  and  two  dimensional  axial  sym¬ 
metrical  two-phase  flow  in  nozzles.  There  has  been  some  progress. 
Two-phase  flow  in  a  combustion  chamber  is  involved  with  mass, 
and  has  some  new  characteristics.  The  author  has  done  research 
on  one  dimensional  two-phase  constant  lag  flow  in  a  combustion 
chamber  [4]  ,  analyzing  the  particle  velocity  lag  effects  on  the 
combustion  chamber  processes.  This  article  studies  the  one 
dimensional  two-phase  nonaquilibrium  flow  in  a  combustion  chamber. 
Based  upon  the  principle  equations  of  one  dimensional  mass  involving 
two-phase  flow  in  a  combustion  chamber. ,  it  discusses  the  numerical 
solutions  of  the  equations  in  detail,  analyzes  particle  size 
effect  on  the  chamber  flow  field  and  the  pressure-time  curve, 
and  compares  the  results  with  that  of  constant  lag  flow. 

Some  practical  conclusions  are  reached.  This  gives  precision 
in  predicting  pressure-time  curve,  and  more  accurate  boundary 
conditions  for  calculation  of  two-phase  flow  in  jet  nozzles. 
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II.  FUNDAMENTAL  EQUATIONS 


Assume: 
1.  Flow 


3. 


is  one-dimensional  and  steady-state; 

2.  The  friction  and  heat  loss  on  the  duct  wall  are  negligible; 

A^O^  particles  are  spherical,  uniform,  and  in  the  liquid 
phase;  the  particle  volume  and  the  Brownian  motion  effects 
on  pressure  are  negligible; 

The  gas  phase  is  an  ideal  gas  at  freezing  point  except  in 
contact  with  particles  and  otherwise  is  inviscid; 

No  mass  exchange  between  two  phases; 

The  specific  heat  of  gas  and  particles  is  constant  /13 


Based  upon  the  above  assumptions,  the  fundamental  equations 
of  two-phase  flow  in  the  combustion  chamber  are  obtained 
as  follows  DJ  and  C3j : 

Gas  phase  - 
Mass  Equation: 

d . 

(1) 


—  e)pTrS 


Momentum  equation: 


-jfctprtAl  —  A^-XprA 


(2) 


Energy  equation: 

d  {  u*  \  1 

Pt9VtAyht  +  A o  y*~)  {  =  (1  —e)pTrS gH sg  —  AoX ppVpA  +  qprA 


dx\ 

in  which 


(3) 


Condensed  phase- 
Mass  Equation: 


r  =  bP" 

h  o  —  C  rtT  $ 


-j^(prVrA)=  epTrS 


(4) 


\  s»» 


Momentum  equation: 


JL.(p,.v''A)-Xp'A 


Energy  equation: 

Af prgvrA(hr+A<)  =  eprrS gH s?  +  AqX pt-vpA  —  qppA 


dX\rra’r"\ - 2  g 

Two-phase  mixture- 
Mass  equation: 


dth  c 

jj-P.rS 


Momentum  equation: 


dP 

jL^p.v'.A  +prv',A)=-A-j^ 


Energy  equation: 


■fi*.  fj)  %)]-‘”rSH: 


in  which 


=ptv,A 
fhP  =p,vFA 
1h  *  fftf  +  thr 


III.  SOME  SUPPLEMENTARY  RELATIONSHIPS 

Equations  ( 1 )  —  ( 6 )  ire  not  closed,  therefore, 
we  introduce  the  following  supplementary  relationships: 


1. 


Gas-phase  condition  equation 

for  an  ideal  gas  —  P=  gR>T 


■.>  *p  w^’^'y,y^M>wM’wwjwww^Tyf¥'Jwyjiwwyjf|jPTy.|w|ji 


2. 


Condensed-phase  condition  equation. 

When  the  Al^O^  particle  temperature  is  greater  than  the  melting 
point  (7\_-23i8*K)  assuming  its  specific  heat  is  constant, 
then 

h?**  h,m  + Ci(T  ?  —  Trm)  (11) 

in  which 

-  Specific  heat  of  liquid  A^O^  particle,  0.34327  KCal/kg. 


degree  K;  h  -  Enthalpy  of  liquid  A^O^  at  Tpm, 


pm' 


876.9498 


KCal/kg 


From  (11)  we  have 

dhp^CtdT,  (12) 

The  particle  resistance  X[l,2]  per  unit  mass 

X  =  A,(v,  —  vt)  (13) 

Under  the  condition  of  combustion  chambers,  the  particles 
carryout  Stokes  flow,  at  this  moment 


Ar  =  -  -&*- 


2  r:,pmr 

2 

in  which  p,  =  l  .208  x  l(T'7'‘/A7’-,(k9/  sec/m  ) 

4.  The  heat  flux  q  [1,2]  of  particles  per  unit  mass 

q-Bp(T,-T9)  (14) 

in  which  B,=  — 

r\p-r 

5.  The  total  enthalpy  Hs  in  1  kq  of  two-phase  mixture  Hs  =  (l  —  e)Hst  +  eHSP 
From  the  assumption  2,  Hg  is  constant  along  the  channel, 

therefore,  it  is  convenient  to  use  the  parameters  0»  =  i>=  0,  Te  =  TP  =  T,) 


3 

V 

V 

s 

X 


at  the 


charge 


origin  (x=0),  to  express  H  .  Therefore 


/7  s  =  ( 1  —  e)C  PtT,  +  e[hPm  +  C  i{T ,  T  >>«)] 


(15) 


a 


V 

3 
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IV.  COMPUTATIONAL  EQUATIONS 


(OsSAT^l) 

Vg 

(0<L<1) 

*  *— 1  » 


Therefore  particle  velocity  lag  = 


Vg  —  Up. 


■  1-K 


particle  temperature  -lag  =  Zj. — Zj.  =*!  —  /, 

T,  —  Tg 

Here  k,L  are  defined  as  particle  velocity  lag  coefficient 
and  temperature  lag  coefficient  respectively. 

Ignoring  the  effect  of  burning  erosion,  assuming  the  cross 
section  area  A  of  the  charge  duct  is  constant  along 
the  longitudinal  axis,  after  an  elaborate  manipulation,  the 
following  numerical  solution  is  obtained  based  upon  the 
fundamental  equations: 

dvP  _a  Vg  —  Vg_.  eprrS  \ 

dx  r  vP  pgA 

dvg _ Rt _ _ _ _ — —  i  gp phgAr-^- — — 

dx  Pg{v',Cpg  —  gRgTgCpg  —  AoRgVt)t-  Vg 

5 

+  pPBp(Tp-Tg)+AQprAp(v,-vpy  +  (\-e)prrjg 


l  u  _l.  A  (1— g>PrrS  _  g  A  p  Vg-Vr 

•  (Hs,+h,+AQ—)\  l-«  ' 


dm  _  c 


v,A 

_(1  -e)m 
Vg  A 

P=*Pt-(p,Vl,+  PrV\) 

j  —  E- 

'  P,gR' 

T P  =  — kr^-  &CiT » +  (1  —  pgT t~ (l- pgT g 

fiC/u 

Vf  £*  - 
JkgRgT, 

T-  T  l ,  1  k—  1  \ 


xn  which 


Hst  =  C„Ts 


The  system  of  equation  (18)  contains  3  ordinary  differential 
equations  and  8  algebraic  equations,  with  the  unknowns. ^  ^ r '  v *' 


P#-.P,  Ps.  MgRth,  un(jer  the  given  boundary  conditions,  which 


can  be  solved  by  the  Runge-Kutta  numerical  method. 


V.  INITIAL  CONDITION  AND  BOUNDARY  CONDITIONS 


The  initial  condition  of  the  equation: 


When  t=0 ,  S=S . ,  A=A.  .  The  value  of  A.,  S-  can  be  obtained  by 
i'  10  xi  1 

the  actual  shape  of  a  charge  duct,  for  a  circular 


cross  section  of  the  duct. 


S  ,  —  2xyt 


y^  is  the  duct  initial  radius. 


The  boundary  conditions  of  equations: 


If  the  charge 


fills  the  annular  space  between  two  closed 


end  coaxial  cylinders,  the  boundary  conditions  are  (refer  to 


Figure  1);  at  the  charge  head  (x=0)  cross  section 


v/’t  =  vl,  =  0 


rh,  =  0 


Ttt  —  Trc  —  T 


P  =  P. 


gR,T, 


1  -e  vr,  1  —  eK, 


I 

! 

in  which,  the  total  temperature  of  the  charge  head  gaseous  phase  Tq  is  equal 
to  the  burning  temperature  of  the  charge..  Pq  is  the  correspond¬ 

ing  pressure.  Before  solving  the  equation,  it  is  unknown.  Therefore, 
its  value  has  to  be  determined  by  the  iteration  process  in  the 
numerical  solution.  The  particle  density  Pf*  at  the  charge 
head  is  ^  an  uncertain  value.  To  determine  this  value,  the 
value  of  Kq  has  to  be  obtained  beforehand. 


At  the  charge  end  (x*l)  cross  section  the  mass  flow  rate  m  of  the 
two-phase  mixture  passing  through  the  cross  section  area  at  the 
charge  end  (x-1),  should  be  equal  to  the  mass 

flow  rate  m  of  the  two-phase  mixture  passing  through  the  throat 
of  the  jet  nozzle,  that  is 

(20) 

ft»L  can  be  determined  by  the  parameters  of 
of  the  x= I  cross  section,  that  is 

(21) 

Ttl  L  =  A\Pq  lV 9  1+  p  r  L^'  p  L) 
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The  present  article  is  mainly  about  two-phase  nonequilibrium 


flow  in  a  combustion  chamber.  To  avoid  the  numerical  solution 
for  one  dimensional  nonequilibrium  two-phase  flow  in  the  jet 
nozzle,  in  the  determination  of  the  boundary  condition 
of  the  charge  end,  it  is  assumed  that  the 

flow  in  the  jet  nozzle  is  a  one  dimensional  two-phase  constant 
lag  flow,  therefore: 

gR,T,C  (22) 


in  which 


C  ~ 1  +  iZ~t K C*(  1  - K)  +  AT]  +  (k- 1) 6LD} 


D  = 


1  + 


1-V 


1+-J-SL 
1  —  e 


1+  (A-l)§ 


«^-i» 


where  the  value  of  K,  L  should  be  taken  as  the  corresponding 
value  at  the  end  of  the  field  length. 


VI.  THE  DETERMINATION  OF  K  ,  L 
(  °  ° 

For  obtaining  the  distribution  of  particle  speed  lag  and 
temperature  lag  along  the  duct,,  and  the  numerical  solution 

for  a  system  of  equations  (18) ,  the  values  of  K,  L  at  the 


charge  head  (x=0)  have  to  be  determined. 

However,  at  x-0 ,  both  K  and  L  are  indeterminate  forms,  therefore, 
it  is  necessary  to  consider  that  at  the  charge  head 

i 

) 

I 

! 

I 

1.  The  determination  of  K 

o 

'  from  1  1  Hopital's  rule 


K. 


(a) 


After  transforming  the  first  equation  in  the  system  of  equations 
(18)  ,  we  find  that  ,  at  x-*0*  \  is  also 

an  indeterminate  form  of  So 


Vp 


'Ap 


v:\~s— vp\^+  _ 

V  r  1  r»l 


For  convenience,  the  subscript  x-»0+  has  been  omitted  in 
the  following  expression.  After  rearrangement  it  reads 


I 


-  v. 


(b) 


With  the  help  of  the  boundary  condition  (19) ,  from  the  second 
equation  in  the  system  of  equations  (18)  the  following  is 
derived  that  at  x-»0+ 


in  which  r.^bPl 


T.>-Q— *W.S 

P..A 


(c) 


Substitute  (c)  into  (b) ,  it  is  obtained  that  at 


x-»0+ 


Because  in  a  combustion  chamber  the  flow  involves 


an  additional  mass  of  particles,  if  the  speed  of  a  particle 
increases,  then  therefore,  (d)  should  be  positive 

which  results  in  the  sign  in  front  of  the  square  root  being 


Substitute  (d) , 


then  obtain 


(c)  into  (a)  ,  and  set 

^  4  (l-t)prr,5 

a  Ar  P»*A 


K.~ 


-1  +  Jl  +  Za 
a 


c 


c 


The  determination  of  L  : 

o 

Same  as  above,  from  1'  Hopital  rule 


It  can  be  derived  from  the  fundamental  equations  ( 1 ) — ( 6 ) 


dTP 

dx 


r-  1  In  —h  4-  A  ^  \  r  —  T,} 

=  +  ij)  VrffC i 


herein,  T'A*—*  is  also  an  indeterminate  form  of 


0 

0. 


From  1' Hopital  rule  it  is  obtained  as 


Tr\s  =  -T’r  | 


B  r  T  f  1 »— t  *  T ,  \  r—,* 
gCi  1 


Omitting  the  subscript  ,  then 


from  equation  (f)  it  is  obtained 


I 


(25) 


With  the  expressions  of  equations  (d)  and  (23) ,  the  following 


is  obtained: 


.  *  B,  1 

P  gC,  Ar  (  —  1  +  %/  1  +  2a) 


(26) 


VII.  NUMERICAL  SOLUTION  OF  THE  SYSTEM  OF  EQUATIONS 

Under  the  boundary  conditions  (19)  ,  (20)  and  the  initial 
condition,  but  use  of  the  fourth  order  Runge-Kutta  method  the 
solution  for  the  system  of  equations  (18)  can  be  obtained  as 
follows : 

1.  At  a  fixed  time  t,  the  gas  flow  parameters  (t>»,  T „  pt,  P, 

Pi,  T j.  Mi),-  ,  the  particle  flow  parameters  Pp)  and 

the  particle  lag  factor  K,  L  distribution  along  the 
duct . 

> 

f 

!  2.  The  pressure  variation  with  respect  to  time  in  the  combustion 

I 

'  chamber. 

In  the  boundary  condition  (19),  the  head  pressure  is  unknown  apriori, 
and  has  to  be  determined  with  the  numerical  solution  process.  The 
procedure  is  the  following: 


(1)  Form  Pi1’— (C*prbK to  compute  the  first  approximate  value; 

(2)  With  the  value  of  P  ^  solve  for  the  system  of  equations 

(18) ,  and  obtain  the  variation  of  parameters  of  the  gas  and  the 
particle  along  the  x-axis,  then  form  equation  (21),  (22)  and 

solve  for  rii^  and  mfc; 

(3)  Let  Artiom,  —  fhL  -,  and  make  a  judgment  on 

|^|<«,  (2’) 

I  iht  I 


(  is  a  given  allowable  error,  for  instance  £^=0.01)  to  see 
whether  it  is  satisfied.  ’ 

(4)  If  equation  (27)  is  not  true,  and  a*»»>0!  ,  take  Pi”— Pi1*—  AP 

(^P  is  the  given  increase  in  pressure,  for  example 

AAP=1  kg/cmx)  ,  otherwise  take  P(,’,=  Pl',+ AP  »  repeat  the  pro¬ 
cedure  for  solving  the  system  of  equations  (18)  ,  if  the  mass 
flow  rate  satisfies  (27),  then  PQ  ^  is  the  true  value  of  the 
pressure  at  the  charge  head  PQ; 

(5)  If  the  condition  (27)  still  can  not  be  satisfied,  the  following 
interpolation  equation  can  be  used  for  computing  P  ^  : 


p< 

*  v 


»<  •). 


pi  >_  p<  •-!  »  „  N 

. c-3'*-3- . > 


then  solve  for  the  system  of  equations  (18)  ,  till  the  mass  flow 
rate  satisfies  equation  (27)  .  Figure  2  shows  the  flow  chart 


of  the  numerical  solution  process,  in  which  'r* 


cp 


burning  speed  along  the  length  of  the  duct. 


is  the  average 


lV- 


Figure  2. 

Key:  (1)  start; (2)  input; (3)  compute  PQ;  (4)  furnish 
boundary  conditions;  (5)  Runge-Kutta  method  for  solving 
the  system  of  equations  (18),  (6)  output,  (7)  adjust 
the  head  pressure  (8)  No,  (9)  Yes, 

(10)  Yes,  (11)  No,  (12)  Stop,  (13)  Flow  Chart  of 
Computation . 
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mw,  -vrwn  »ngr.  wk\v,  •rjw.nr. 

VIII.  THE  EFFECT  OF  THE  PARTICLE  SIZE  ON 
INTERNAL  BALLISTIC  PROPERTIES 

As  to  a  certain  solid  propellant  rocket  motor  (  e=0.26), 
the  pressure-time  curve  and  the  flow  field  in  the  combustion 
chamber  have  been  computed  for  different  sizes  of  particles. 

Figure  3  shows  the  effect  of  the  particle  size  on  the  Pressure- 
Time  curve.  It  can  be  observed  in  the  figure  that  the  pressure 


in  the  combustion  chamber  decreases  as  the  particle  size  increases. 

At  t=0  sec. ,  the  variation  of  the  combustion  chamber  pressure 

PQ  with  respect  to  the  particle  radius  is  shown  in  Figure  4. 

The  decrease  of  the  combustion  chamber  pressure  will  lead  to 

a  slowdown  of  the  burning  speed,  increase  the  charge  burning 

period  and  decrease  the  mass  flow  rate. 

At  t=0  sec.,  the  variation  of  the  gas  parameters  along  the 

x-axis  for  the  different  particle  sizes  is  shown  in  Figure  5. 

Mg  solely  increases  with  the  distance,  the  rest  of  the  parameters 

decrease  with  the  distance.  Because  PQ  decreases  when  the  particle 

size  increases,  therefore  P,P  ,  p  are  affected  by  the  particle 

s  g 

size  distinctively.  The  larger  the  particle  size,  the  smaller 
are  the  values  of  them.  But  the  particle  size  effect  on  Tg  is 
negligible . 

Figure  6  shows  at  t=0  sec.  the  distribution  of  the  gas  speed 
,  g  and  the  particle  speed  Vp  along  the  duct  under 
different  particle  size  conditions.  Because  the  flow 
in  the  combustion  chamber  is  involved  with  an  increase  in  mass, 
both  v-g  and  increase  along  the  length  of  the  duct.  But  the  particle 
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size  effect  on  \r  and  \f  is  different  when  particle  size  increases 

g  p 

\Tg  increases,  but  decreases.  This  is  because  the  particle 
of  unit  mass  has  a  resistance  to  the  gas  of  92  p-  .  When  the 

'  P. 

particle  size  increases,  X  decreases,  so  the  gas  speed  increases. 

At  this  moment  the  torque  exerted  on  the  particle  by  gas  decreases. 
This  results  in  the  particle  acceleration  decreasing,  therefore 
'O'  becomes  small.  The  density  and  temperature  distribution  of 
the  particle  are  shown  in  Figure  7.  It  can  be  observed  that 
£p  hardly  decreases  along  the  duct  and  T^  basically  stays 

constant. 

The  effect  of  the  particle  size  on  the  value  of  is 
significant:  when  particle  size  increases,  increases  too. 

But  the  effect  on  the  value  of  T^  is  negligible. 

At  t=0  sec.,  the  particle  speed  lag  factor  K  and  temperature 
lag  factor  L  distributions  along  the  duct  are  shown 
in  Figure  8.  It  can  be  observed  that  the  value  of  K  slightly 
decreases  along  the  distance,  but  the  value  of  L  fundamentally 
stays  constant.  This  indicates  that  in  the  combustion  chamber, 


the  particle  speed  lag  (1-K)  increases  slightly  along  the  length  of 
the  duct,  but  the  temperature  lag  (1-L)  basically  stays  constant. 

The  smaller  is  the  particle  size,  the  smaller  is  the  variation 
of  K  value  along  the  distance.  (Refer  to  Table  1.)  Therefore 
when  the  particle  size  is  small,  the  two-phase  flow  in  the 
combustion  chamber  can  be  treated  as  the  constant  lag  flow  motion. 


L’TWW^  VWTfVTjvVl ITT  ^ Vw  v*  yr* v  r  x *  w  * W “ *  V  'P  'V  r:w~. 
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Table  1. 


Key:  ( //  The  radius  of  the  particle  r 

It 

(micrometer) ,  (  »  The  decrease  of  K  at  the  charge  end 

Table  2  lists,  for  different  particle  sizes,  the  comparison  between 
the  numerical  solutions  and  constant  lag  flow  computations  of 
some  important  parameters  in  two-phase  flow  in  a  combustion 
chamber.  Compared  with  the  numerical  solution,  the  combustion 
chamber  pressure  PQ  and  the  total  pressure  at  the  end 
charge  end  obtained  from  the  constant  lag  flow 

computation  are  comparatively  low,  and  the  rest  of  the  parameters 
are  comparatively  high,  but  the  deviation  is  less  than  1%. 


Table  2. 


2  2 

Key:  (2)  Micrometer,  (3)  Kg/cm  ,  (4)  kg/cm  , 

(5)  m/sec,  (6)  m/sec,  (7 )  ( 9),(  11 ),( 13)  numerical  solution, 

(8), (10). (12), (14)  constant  lag  flow. 


From  equation  (24) ,  (25)  it  becomes  known  the  values  of  KQ,  LQ 

are  functions  of  the  particle  radius  r  .  Their  variation  with 

respect  tor  is  shown  in  Figure  9.  It  can  be  seen  in  the  Figure 
P 

K  ,  L  values  decrease  when  the  particle  size  increases,  which 
o  o 

means,  the  particle  speed  lag  and  temperature  lag  increases  as 
the  particle  size  increases. 

Besides,  Kq  <  Lq,  in  the  nonequilibrium  flow  in  the  combustion 
chamber,  the  lag  of  the  particle  speed  is  greater  than  the  lag 
of  the  particle  temperature. 

IX.  CONCLUSION 

Through  the  above  discussion  the  following  conclusions  are 
therefore  obtained: 

1.  Because  of  the  effect  of  the  two-phase  flow,  the  pressure 

of  the  combustion  chamber  decreases,  and  the  larger  the  particle 
size,  the  smaller  is  the  combustion  chamber  pressure. 

2.  The  two-phase  flow  affects  the  flow  field  in  the  combustion 
chamber  greatly.  When  the  particle  size  increases,  the  gas 
speed  increases,  the  pressure,  total  pressure  and  the  gas 
density  all  decrease.  When  the  particle  speed  decreases, 
the  density  increases.  The  effect  of  the  particle  size  on 
the  temperature  of  both  gas  and  particle  is  very 

small . 

3.  In  the  combustion  chamber  the  lag  of  the  particle  speed  increases 
along  the  length  of  the  duct,  the  lag  of  the  temperature  basically 
keeps  constant.  Besides,  the  particle  speed  lag  is  greater 
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SUMMARY 


On  the  basis  of  the  evidence  presented,  a  steady  state  combustion 
model  (l)  of  solid  state  propellents  with  negative  pressure  coefficient 
burn  speed  characteristics  allows  us  to  derive  a  new  pressure  response 
function  formula.  This  can  be  used  to  explain  pressure  pairing  phenomena 
for  propellents  whose  combustion  speed  pressure  exponents  are  zero, 
positive,  and  negative.  The  burning  propellent  is  divided  into  two 
sections:  the  first  section  is  a  structure  composed  of  oxidizing,  agent 

covered  by  nolten  binding  agent  and  the  corresponding  binding  agent 
^i-face.  The  other  section  is  formed  by  uncovered  oxidizer  surface  and 
remaining  binding  agent  surface.  This  model  is  different  from  the 
various  types  of  models  in  the  past.  In  the  combustion  on  the  surface  of 
the  first  type  of  section  described  above,  consideration  has  been  given 
to  the  oxidizer,  under  conditions  in  which  it  is  covered  by  molten 
binding  agent,  so  that  it  is  considered  to  be  in  a  state  of  opposed 
gasification  and  congealed  phase  reaction.  Therefore,  the  real  section 
of  the  pressure  response  function  which  is  obtained,  when  the  pressure 
exponent  of  the  propellent  steady  state  fuel  speed  is  zero  or  has  a 
negative  value,  is  also  capable  of  being  a  positive  value.  When  we  made 
use  of  the  expression  obtained  for  the  pressure  response  function  in 
experimentation  with  the  propellent  (S04-5A)and  made  qualitative 
calculations,  the  results  of  these  calculations  satisfactorily  explain 
the  phenomenon  of  propellents  with  negative  pressure  exponents  still  being 
unstable  in  co  abas t ion  when  aoet  of  the  surface  area  of  the  oxidizer  is 
covered  by  molten  binder  agent.  This  not  only  overcomes  the  weaknesses 
which  all  expressions  for  pressure  response  functions  had  in  the  past 
when  used  with  negative  pressure  exponent  propellents,  but  also,  in  a 
way,  reflects  the  accuracy  of  the  solid  fuel  propellent  steady  state 
combustion  moiel  (l)  in  its  combustion  speed  characteristics  for 
positive  and  negative  pressure  exponents. 

Explanation  of  Symbols 


A 


4 
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■V.V--A 


Oxidizer  Gasification  Rate  Indicator  Prefactor 


A 


•* 


B 


_  2R*T*i 


Indicator  Prefactor  in  Oxidation  Agent  Equilibrium  Evaporation 
Te  Formula 


Pexp  (*%) 


c  opecific  Heat 


E 


C+q/E\m 


O+A  t(AB-A) 


E,.  Energy  of  Activation  for  oxidizer  Agent  Gasif ication  Reaction 

E/  Energy  of  Activation  for  the  Propellent  Gas  Phase  Combustion 

•)  o  s  s  a 

p  bA+  A,(AB  —  A) 


G  Hass  fraction  of  Oxidizer  Agent  for  Completion  of  Condensation 

i  -  -  -'eaotion 

K  ".peal  Constant  for  Gas  Phase  Reaction 


m  Mass  Plow  Rate 

„  Combustion  Speed  Pressure  Exponent 
p  Pressure 

Qf  Unit  Mass  Binder  Decomposition  Reaction  Heat 
0.  Boundary  Layer  Reaction  Heat  of  Unit  Mass  Oxidizer 
9  Heat  of  Evaporation  for  Unit  Mass  Oxidizer 

<h  Heat  Plow  for  Gas  Phase  Boundary  Surface  Plow 
0»  Amount  of  Heat  Released  by  Unit  Mass  Gas  Phase  Reaction 

q<  Rate  of  Release  for  the  Heat  of  Purification  from  Boundary  Surface 

Reaction 

0.  Real  Section 

/?'  Universal  Gas  Constant 


Rr  Pressure  Response  Punction 


ml  P 


S  =iaii_=ffl  (  i  is  the  imaginary  number  unit  J— 1  ) 

m'c 

j  Temperature 

Tt  Absolute  Temperature  of  the  Propellent  Flame 
Ti  Initial  Propellent  Temperature 


T,  °r*opellent  Surface  Temperature 


Wam>  rram  molecular  weight  of  oxidizer  evaporation  gas 

WG  Gram  Molecular  Weight  of  Gaseous  Products  from  Congealed  Phase 
Oxidizer  Reaction 


x  Distance 

a  -lass  “Fraction  of  Oxidizer  in  the  Propellent 

y  Surface  Area  “Fraction  Composed  of  Oxidizer  Agent  Surface  Covered  by 
Molten  Binder  Agent 

X  Congealed  and  Solid  Phase  Heat  Conduction  Coefficient 
i  Hondimensional  Distance 


p  Congealed  and  Solid  Density 


G  Hoad imensional  Frequency 

IB  C 


co  Angular  Frequency 


SUPERSCRIPTS 


—  Steady  State  Value  or  Average  Value 


,  ,Js"lu”bet  “on  Value 


~  Co  ial  ->  <  \  ar>l :  tud  e  of  Perturbation 


.'■■.y.v.v a-'/S-y. 
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1  Parameters  Corresponding  to  Area  I 


1  Parameters  Cor  responding  to  Area  II 


3U33CRI PTS 


«  imaginary  Section 
r  Real  Section 

I  Parameters  Corresponding  to  Area  I 

1  Parameters  Corresponding  to  Area  II 

s  Soundary  Surface,  Value  (Congealed  Phase  Side)  for  the  Place  Where 
*-o 

1  ,  Introduction 

Sonic  instability  is  the  primary  form  of  combustion  instability  of 
solid  propellent  rocket  motors,  ’doreover,  this  type  of  sonic  instability 
has  as  its  primary  source  the  combustion  response  of  solid  fuel 
prop -11  ants .  Because  of  this,  in  the  test  construction  of  solid  rocket 
actors,  the  sonic  instability  which  they  show  is  extremely  important,  and 
the  theory  and  measurements  associated  with  combustion  response  are  I 

indispensable.  In  linear  sonic  analysis,  on  sees  shown  whether  or  not  a  j 

soli  1  -ooket  motor  will  give  a  timely,  consistent  pressure  perturbation  j 

a  a!  irg  -  -i  eot ,  forming  a  sonically  unstable  fuel.  In  this  type  of 
;  i  tty  ’  i  ■ ,  the  pressure  r  espouse  function  must  already  be  known.  1 

"1--'  ~  “n  1  o  '  )bl 3  l  o p  solid  propellent  pressure  response  has  received 

l  i  'gj  a  r.o  jint  s  of  research  work  in  both  the  U.3.  and  the  Soviet  Union. 

! 


F.E.C.  Culick  (2)  and  N.S.  Cohen  (3)  have  already  done  excellent 
critiques  of  this  problem.  Concerning  "flatbed”  propellent  pressure 
response  functions,  F.E.C.  Culick  and  others  (2)  (4)  have  also  done  large 
amounts  of  research.  However,  the  author  recognizes  that  the  steady 
state  combustion  model  and  several  assumptions  which  form  the  basis  of 
the  research  are  still  short  on  experimental  data  to  support  them. 
Moreover,  several  points  of  theory  are  still  doubtful.  For  example,  in 
order  to  explain  the  phenomenon  of  sonic  vibrations  in  the  combustion  of 
flatbed  propellents,  we  take  the  mass  gasification  rate  for  solid  phase 
reaction  areas  and  change  it  for  use  in  the  expression  below: 

f. 

m  =  A,P‘>e 

From  this,  we  can  derive  the  pressure  response  function  formula: 

r  __  nAB  +  n,(A—]) 

A  +  ^-(l  +A)  +  AB 


Moreover,  when  we  use  this  relationship  to  compare  with  the  results  of 
■i  xper  I  non  ts ,  we  select  n,  =  i.o  .  However,  we  know  that,  when  the 

pr  i  i  Mrs  exponent  of  the  auss  gasification  rate  for  solid  phase  reaction 
areas  is  a  positive  value,  the  pressure  exponent  for  its  steady  state 
combustion  speed  can  hardly  be  zero.  This  article  does  not  intend  to 
make  more  criticism  of  these  publications .  We  will  present  on  the  basis 
of  a  steady  state  model  (1)  of  solid  fuel  propellents  in  terms  of 
positive  and  negative  pressure  exponents  for  conbustion  speed 
character i st  ics ,  and  derive  a  new  pressure  response  function  formula  and 
fo-’oe  diagram.  We  will  use  the  i  to  make  a  complete  explanation  of  the 
pr  as  sure  pairing  pheno  nenou  w’n'c’n  exists  in  various  types  of  solid 
propellents  when  their  combustion  speed  pressure  exponents  are  zero, 
positive,  and  negative. 


2,  Physical-Chemical  Model 

In  the  same  way  as  was  the  case  with  the  corresponding  steady  state 
model  (l),  we  select  for  use  a  one  dimensional  aodel  and  divide  the 
combustion  surface  of  the  propellent  into  two  sections.  One  section  is 
exposed  of  oxidizer  agent  surface  covered  by  molten  binder  agent  and  the 
o  v'-'e -sp  )nl  lag  binder'  agent  surface.  Let  us  also  stipulate  that  the 
fraction  of  surface  area  which  it  occupies  out  of  the  whole  combustion 
surface  is  v  .  The  other  section  is  then  composed  of  oxidizer  surface 
which  is  not  covered  and  the  remaining  binder  agent  surface.  These  two 
sections  are  respectively  called,  for  short,  Area  I  and  Area  II,  as  shown 
in  7lg.  1  .  At  the  same  time,  we  take  the  combustion  process  for  each 
solid  propellent  area  and  simplify  it  into  three  stages,  which 
respectively  occur  in  three  different  phases  (Pig.  2). 

1)  The  solid  phase  area  on  the  inside  of  the  solid  propellent,  to 
which  heat  is  applied,  is  gasified  by  the  gas  phase  flame  area  with  the 
congealed  phase  reaction  area  supplying  the  heat. 

2)  Within  the  congealed  phase  layer  between  the  solid  phase  and  the 
gas  phase  one  sees  the  development  of  an  oxidizer  boundary  layer  reaction 
and  binder  agent  heat  of  decomposition  which  are  contained  in  the 
congealed  phase  reaction  and  the  gasification  reaction.  In  Area  I, 
unlike  Area  II,  one  must  consider  the  molten  binder  agent  liquid  layer 
which  covers  the  surface  of  the  oxidizer  as  well  as  the  oxidizer 
underneath  It  as  they  exist  in  the  congealed  phase  reaction  related  to 
gasification  and  its  opposite  process — opposed  gasification. 

3)  In  the  gas  phase  one  sees  the  occurrence  of  combustion  processes 
which  include  dispersion,  mixing,  and  chemical  reaction. 
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Pig.  1  One  Dimensional  Steady  State  Model 
1.  Dispensed  Plame  2.  Perturbed  Flame  3*  Oxydizer  Agent  4.  Liquid 
Phase  5-  Solid  Phase  6.  Paired  Agent  Mixed  Layer  7.  End  Plame  8. 

AP  Plame  9.  Initial  Plame  10.  Combustion  Surface  11.  Propellant  12 
Liquid  Phase  13*  Solid  Phase 


/ 
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Pig.  2  4  Simplified  Two  Dimensional  Model  1 .  Congealed  Phase  Layer 

Solid  Phase  4rea  of  Increased  Heat  3*  Sas  Phase  Reaction  4rea 


3,  Mathematical  Treatment 


For  the  sake  of  convenience  in  the  mathematical  treatment,  we  first 
make  the  following  assumptions: 

1  )  The  influence  of  heat  damage  is  negligible. 

2)  The  solid  propellent  is  noncompressible ,  of  a  uniform  qualtiy, 
and  of  the  same  nature  in  all  directions. 

3)  The  congealed  phase  layer  between  the  solid  phase  and  the  gas 
phase  is  an  infinitely  thin  plana.  In  this  plane  one  sees  securing  the 
oxidizer  congealed  phase  reaction  and  the  gasification  reaction  (Area  I 
and  Area  II  are  different  and  each  has  its  own  rules  of  reaction) .  At 
the  same  time,  one  sees  occuring  the  high  temperature  decomposition  of 
binder  agent.  This  plane  is  called  the  boundary  surface. 

(4)  In  the  three  solid,  liquid,  and  gas  phases  discussed  above,  the 
specific  heat  is  always  a  constant. 

5)  las  phase  reactions  can  be  seen  to  be  steady  stats  processes. 

6)  The  influence  of  congealed  phase  reactions  on  pressure 
disturbances  is  negligible. 

7)  When  the  amount  of  perturbation  is  sufficiently  small,  its  second 
and  higher  order  small  amounts  can  all  be  neglected. 

At  the  same  time,  in  the  same  fashion  as  normal  models,  we  choose 
for  use  a  moving  coordinate  system.  These  coordinates  take  the 
instantaneous  burn  speeds  of  the  propellents  and  move  them  into  the  solid 
phase  area.  In  this  way,  the  origin  of  the  coordinate  system  cun  fall 
entirely  on  the  solid-gas  phase  boundary  surface. 

Below  we  present  the  respective  treatments  for  Area  I  and  Area  II. 

Area  I : 

(l)  Basic  Equations 

1 )  The  unstable  thermal  conductivity  equation  in  the  solid  phase 
heat  addition  area  is: 


.  d'-T  17’  ?>T 

AW‘~mc  aT=pc  * 
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2)  When  we  consider  the  existence  of  molten  boundary  agent  on  the 
cohering  of  oxidizer  surface  and  the  oxidizer  under  the  covering  as  they 
affect  congealed  phase  reactions  and  gasification  reactions,  from 
reference  article  (l),  the  mass  gasification  rate  of  propellent  is: 


3)  If  ye  assume  a  steady  state,  the  mass  combustion  rate  for  the  gas 
phase  reaction  area  can  be  expressed  as  (5): 


i-KP"ex p(  2R*T  / ) 


4)  Energy  equilibrium  equation  for  the  inside  of  the  infinitely  thin 
congealed  phase  layer: 

The  difference  between  the  heat  flow  <?«  from  the  gas  phase 

,dT  1 

combustion  area  and  the  heat  flow  A!>x\  transfered  into  tne 

solid  phase  should  equal  the  rate  of  heat  release  qr ,  which  is  produced 
by  the  reaction  in  th?  congealed  phase  layer  of  the  congealed  phase  and 
the  gasification  reaction,  that  is, 


Moreover,  ;  *  l  ug  motion  coo  rd  i  nates  for  explanation,  Qi  =  mLQt-c(T/—T,)2, 

,  also,  qc  =  m[aQ,-(l -a)Qf]o  .  In  these  equations,  Qc  is  the 

amount  of  heat  relent  1  by  the  unit  mass  gas  phase  reaction.  is 

the  algebraic  sura  of  the  effective  heats  of  the  oxidizer  agent  congealed 
phase  reaction  and  the  gasification  reaction.  Qr  is  the  effective  heat 
"rom  the  hi  gh  temperature  thermal  decomposition  of  binder  ugint.  I  ?  ■:  i 


substitute  the  equations  above,  we  get 


+m[aQ.-(\-a)QF-Q,+c(T,-T,)]  =0 

OX 


)  Linearization  Treatment: 

1)  Thirst,  we  carry  out  a  linearization  treatment  of  equation  ( 1  )  . 
We  postulate: 

ra(/)=*l  +  i?,{me,“')  : 


7(*.  0  =  T(x)  +  tf.m*Vf*"}  ' 


If  w e  substitute  in  equation  (l),  and  eliminate  the  high  order  small 
quantities,  as  well  as  making  use  of  the  steady  state  condition 


,  d'T  _  dT  _  0 

a  —  —  #nc  -j - « 

dx*  dx 


n,  e  pint  ion  .1)  changes  to  become: 


Rj^e^'^-N/'me'-'c^-K.^me^c^^RApcTxcoe--}, 


•i  i-;e  r.-i  fn-  as  -tr  *  set  up  for  all  tines  1  ,  therefore,  one  nas 


(7) 


dx'  dx  dx 


tm!f. 
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Then  equations  (6)  and  (7)  can  respectively  be  written 


d^__dT^0 

dv  di 

d'T  _df  «  dT_ 

dV  di  m  di 


2)  In  the  same  way,  in  the  case  of  equation  (2),  one  can  set 


m  =  +  =  +  R,{ihe<m' )  i 

T.-T.  +  T\-T.  +  R.{T*,m>) 

p  =  P  +  P'  =  P  +  R'{Pe‘m'} 


Moreover,  one  can  use  the  steady  state  condition 

ii,/  £..j 

m=  c  ’  1^77  i4»"exp  \  i?°T.  /  1  n  /.  ,  g  W^ijmA 

"'V1  +  l-G  W c  ) 

Division  with  equation  (2)  gives  one 


(10) 


(11  ) 


WWVIV 


1  +  —  =exp yR.Tt  •  Tt  ) 


In  the  equation 


/  G  W 

wc  )_i 

Peycp(-^r) 


E  T'  t' 

Also,  bee '.use  of  the  fact  that  piijr  •  -jr  and  ~Ff  '  f~  are  both 

very  small  (tP?~  ari<^  are  both  the  orber  of  magnitude 

of  10).  Moreover,  is  then  in  the  range  KT'-lO"'  •  Therefore 

one  has 


t: 

'  T, , 

)  FT. 

t: 

*  f. 

“K  r-t; 

JV  \ 

T7, 

)  ~ 1  _  ~fT, 

t: 

'  T. 

Substituting  in  e  q  ration  (12)  and  eliminating  the  smaller  quantities  of 
higher  order,  one,  then,  obtains 

,  „  p' _ q_  r 

T.  E*>  Ll  ?  P*T.  T. 
m  FT,  '  T,  c 


As  shove,  because  of  the  fact  that  all  of  the  equivalent  forms  for  time.' 

•ho  ill  o  •  ••  e t  io,  one,  eonse  q.iently  has 


f  r 


C+TT-  T.  FT.  Pi 


3)  After  one  carries  oat  the  same  type  of  linearization  treatment 
on  equation  (3),  one  obtains: 


A  P 
Tj_  A~"P 

r,  TT 

2  R*T, 


(14) 


4)  In  the  same  way,  equation  (4)  can  be  changed  to  become; 


(15) 


(3)  Solution  of  Equations 

Because  equation  (9)  and  its  boundary  condition  are: 


d'T  dT 
dV  di  =° 


t=— OO  time7  =  7', 
'£*=0  time  f  , 


it  follows  that  its  solution  is  obviously  T**T,+(.T,—T,)et 

dT 

Moreover,  ^±rx=(T,-T,)et  .  If  we  take  this  expression  and 

di 

substitute  it  into  equation  (15),  we  then  get 


(16) 


e  substitute  equation  (13)  and  equation  (10),  equation  (10) 


changes  into 


djT  dT  £../  T,\  r/C+E.,\L-  _£L  £  r< 

dv  ~di  ST~  R'\l  T,  I  LV  C  Jr.  C£„  P/ 


Its  boundary  conditions  are: 


£=  —  cotime  T  =  0 
time  ^  =  T , 


Therefore  a  general  solution  for  the  equation  is: 


7 


c'e 


T,  FT. 
T.  CEt, 


j  _  £*.  7\\  , 

In  the  equations  /?*T,\  7,  /  •  Moreover,  is  then 

determined  from  the  characteristic  equation  x'-X~S  =  0  .  Because, 

when  s  ~  oo  f  T  =  0  ;  therefore,  one  only  has  x -(1  +  •Si  +  iS  )/2 

as  a  characteristic  limit.  Again,  from  the  boundary  condition  £=o 
,  when  T  =  t,,  ,  it  is  possible  to  obtain  a  solution  to  the  equation  as 

follows : 


T 


T.  R*T, 
T.  CE M 


(eU_et)  1 

-J 


Substituting  in  equation  (16),  one  then  obtains; 


,  At  c+-e.._T^  FT.  P\ 
*+S\— C-  T.  CE„  P)  ^ 


(5h 


f  1  71 


Prom  equation  (13),  one  has 


c+e7. 


c-+i- 

ffi  P 


(18) 


If  one  takes  equations  (14)  and  (18)  and  substitutes  them  with  equation 
(17),  one  then  obtains: 


C+- 


’  •*_ 

P 


C~+  „ 

m  P 


r,  <7  P 

c+e:.  ( m~ny  fj_ 

m  P  2R*f, 


C+E„ 

C^  +  ? 
R*T, 
E.. 


1 


2R*T) 

If  we  make  ~  (T.-T^E/  «  ’and  we  make'  use  of  the  pressure  response 

function  definition  Rr= ,  then,  the  above  equation  can  be 

arranged  so  that: 


nAB— 


A  —  1 


C  +  q/E„ 


T - A  +  AB+  — —  .  r-  - 

/  C  +  q/E  „ 


U-l) 


Area  II: 

Based  on  a  method  similar  to  the  one  above,  we  only  take  basic 
equation  (2)  and  change  it  to  be  m  =  A„ex  p(--j~*±-)  •  It  is  then 

'.ble  ho  )b  t  a  i  n  an  expression  for  the  ores  nice  response  function  for 
Area  1 1 ( 5 )  : 


/?n=  — - - - 

A-A  +  AB+  (A-l) 


Again,  on  the  basis- of  the  normal  assumption  that  there  is  a  direc 
proportion  between  combustion  surface  gain  and  combustion  surface  area, 
the  pressure  response  function  for  the  entire  combustion  surface  formed 
by  Area  I  and  Area  II  should  be: 


R,=y- 


X—  1 

n<A>B>--c+T[EZ, 

C 


-f-A,+A,B,+  c+g/E,. 


(A-l) 


+  (l-y) 


n,iAjiBi 


A'  -An  +  A.'Bu  +  CX-l) 


Actually,  the  imaginary  section  is: 


?.,n ,  A  i  B i  D + (x(n  i  A ,  B i  —  Q  +  q/E^ 

(tf  A  =  V  W+E*  * 

,  n  xnM^iiBn(A,D  +  X,F) 

+  u  Y)  D'+F 1 

—  Xrn,Ai  B,E+ (x,n,  A i  B i  — C  +  q/E^  ^ 
(Pr)t=*Y  D'  +  E'  5 


D  =  A  +  lr(AB- A),  E  s 


^_p+x^ab-a). 

F  -^Q  +  X,(AB-A)i 

A'-i!,+Vl!"i+(i+i6fl’)i]ii 

a'-W1[-,+(,+i60,)^ 


4  Results  of  Calculations  and  Discussion 

Waking  use  of  the  steady  state  model(l),  the  various  results  of 
intitial  calculations  on  S04-5A  experiments  (Table  l)  form  the  Initial, 
data  to  substitute  into  the  formula  for  the  real  section  of  the  pressure 
response  function,  yielding  results  such  as  those  shown  in  Pig.  3- 


f\g.  3  7’ne  Real  Portion  of  Two  Low  Pressure  Response  Functions  for 
304-5  4 


The  formulas  obtained  -md  the  results  of  calculations  all  clearly 
show  that  the  real  sections  of  pressure  response  functions  for  "plateau" 
and  "mesa"  propellants  whose  pressure  exponents  are  aero  or  negative  are 
all  capable  of  being  positive  values  larger  than  zero  within  quite  a 
large  range  of  frequencies.  From  the  standpoint  of  theory,  this  explains 
very  veil  facts  associate!  4 i th  severe  phenomena  of  unstable  combustion 
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and  the  existence  of  these  two  types  of  propellents(4) (6^(3) .  Moreover, 
just  as  this  author  calculated  in  reference  (1),  due  to  the  fact  that  the 
abnormal  combustion  of  oxidizer  covered  with  molten  binder  agent  takes 
the  place  of  partial  shut  down,  in  this  way,  under  the  pressure 
associated  with  the  existence  of  widespread  covering,  it  is,  on  the 


contrary,  e 

oscillation 

as y  to 

.  It 

see  the  appearance 

is  then  possible  to 

of  the  results  of  s 

give  a  reasonable 

elf-e  xcite 

explanation 

*  1st 

<r. 

■  ft  ffi 

:  ,  i  « 

a  a  <2 

*  pfti 

;  JTffl  2 

1 

B,  ; 

kg/cm’ 

3.  48  X  10* 

? 

kg/cm’ 

23 

1  35 

c 

Cal/g  c  K 

0.3 

]  G 

i 

0.3 

0.2 

£«. 

Cal  mol* 

2.2x10*  ; 

V 

j 

0.6 

0.8 

E, 

Cal  mol* 

;  6.0X10* 

1  " 

0.572 

-0  693 

i 

fl  ; 

Cal.  ir.ole 

2.08x10* 

!  TSI 

•K 

1038.6 

j  1098.8 

r, 

"  K 

i  2900 

"/ 

0.75 

i  -0.783 

T, 

‘K 

293  jj  Ttn 

•K 

1086.6 

;  iioB.o 

i 

W  A,ttt  g  mol*  117.5  1.  n  / 1 

0.486  1  0.424 

HI  g'mol*  1  28.4  1 

i 

'  Cal, 'cm  sic  °K  1  0.3x10"'  j| 

| 

r  g  cm'  1.665 

1 

Table  1  1  .  Parameter  2.  Unit  3-  Numerical  Value  4.  Sample 

T  simulation 

Vhat  the  dotted  line  in  Pig.  3  shows  is  that,  in  reference  ( 1  ) ,  S04-5A 
propellant,  with  a  large  surface  area  covered  and  a  pressure  of  /J  =  23kg/  cm’ 
,  in  experiments  with  the  "T"  type  engine,  shows  the  appearance  of 
self- excited  oscillation.  The  frequencies  involved  correspond  to  the 
real  sect  ion  values  of  the  r  ■>  spouse  functions  for  q  1 1 1  e  high  pressures. 

^  -»ca  is  -  o "  this,  the  appearance  of  self-excited  oscillation  in  these 
experiments  cannot  be  considered  surprising. 
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Most  of  the  results  discussed  above  were  obtained  through  analysis 
in  4rea  I.  This  was  done  through  the  use  of  our  steady  state  model  and 
with  consideration  given  to  the  covering  of  the  oxidizer  surface  with 
molten  binder  agent  and  to  the  existence  of  a  congealed  phase  reaction  as 
well  as  opposed  gasification.  Because  of  this,  the  theoretical  formula 
derived  for  the  pressure  response  function  can  not  only  overcome  the 
deficiency  associated  with  the  fact  that  the  pressure  response  function 
formulas  we  already  have  cannot  be  used  with  negative  propellent  pressure 
exponents,  but  also,  in  one  respect,  it  reflects  the  vitality  of  the 
stable  state  combustion  model(l)  for  solid  state  propellents,  which 
presents  the  combustion  speed  characteristics  for  positive  and  negative 


pressure  exponents. 


The  formula  which  is  derived  when  the  covered  surface  area  fraction 


is  zero,  then  becomes  the  "two  parameter  formula"  for 
propellents  with  binder  agents  which  are  difficult  to  melt: 


A  +  —  —  ( 1  +  A)-h  AB 


Uni  ig  .is  ^  of  the  unified  expression  we  already  have  for  the  "gas  phase 
quasi-steady  state,  uniform  propellent  unidimensional  model" ( 2 ) ,  this 
article  consequently  derives  a  formula  which  is  capable  of  having  even 
broader  applicability. 


This  article,  before  deriving  this  formula,  made  no  small  number  of 
empirically  quite  significant  assumptions.  for  example,  the  assumption 
that  the  congealed  phase  layer  which  gives  rise  to  the  congealed  phase 
^action  and  the  gasification  reaction  is  an  infinitely  thin  plane  is 
certainly  not  a  good  approximation^).  The  assumption  that  tne  gas  phase 
a 3  3ua  ?s  a  steady  state  has  an  error  less  than  i 0%  only  when  the 


v;  ^ 


j 
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f re^uencies  are  smaller  than  10,000  Hertz(2)(5).  The  compressiblity  of 
propellents  also  has  an  influence  on  response  functions  which  is  within 
10#(9)*  Heat  radiation  losses,  under  certain  conditions,  also  have  an 
obvious  influence  on  response  functions( 1 0) ( 1 1 ) .  As  far  as  congealed 
phase  reactions  are  concerned,  in  our  stable  state  model,  they  are 
basically  pressure  functions( 1 ) .  However,  we  also  conveniently  saw  from 
our  mathematical  treatment,  that,  when  their  influence  on  pressure 
diminishes  ,  error  is,  naturally,  even  more  unavoidable.  Although, 
except  for  congealed  phase  reaction  problems,  the  majority  of  assumptions 
are  obvious,  they  are  also  universal  in  their  application.  However,  this 
explains  the  fact  that  the  formulas  derived  above  are  still  awaiting 
refinement . 

At  the  same  time,  formulas  derived  on  the  assumption  of  homogeneity 
can,  of  course,  only  be  used  in  situations  where  one  has  homogeneity. 
However,  due  to  the  fact  that  propellant  impurities  have  already  been  (57 

used  by  several  scholars  (12)  in  methods  such  as  the  one  that  follows, 
that  is,  taking  a  multi-mode  composite  propellant  and  viewing  it  as  a 
finely  dispersed  oxydizer  granule/fuel  "surface  match  up"  of  random 
arrangement.  When  one  assumes  that  each  "pairing"  is  mutually 
indpendent,  then  it  also  b scones  possible  to  take  the  propellant  surface 
and  rearrange  it  into  a  family  of  hypothetically  dispersed  units  of 
propellant-  This  is  what  is  called  a  hypothetical  propellant.  And,  the 
pressure  responses  of  these  hypothetical  fuels  are  calculated  using 
formulas  deduced  from  homogeneous  theory.  Moreover,  from  this  one 
obtains  the  pressure  response  for  the  whole  composite  propellant.  Just 
as  is  pointed  out  in  reference  (3),  although  this  method  has  the 
drawbacks  of  two  different  models  in  its  application,  it  is  still 
currently  being  followed.  Because  of  this,  this  model  certainly  does  not 
fail  to  explain  pressure  pairing  problem  values  for  composite  solid 
propellents  because  of  its  adoption  of  a  homogeneity  hypothesis. 


lonclusions 


ij  i;  f  )  ^  t  ?  S  Judy  state  ooabustion  uodel  for  positive  and 


negative  pressure  exponent  combustion  speed  characteristics  for  solid 
propellants,  we  deduced  a  new  pressure  response  function  formila.  This 
^ormula  explains,  in  a  reasonable  way,  phenomena  associated  with  unstable 
•nimbus  tion  which  exists  in  all  solid  propellants  with  positive,  negative 


zero  pressure  exponents. 


A  .V  ■-  S 


44 


a-.  IP  mJmJ>  "»  M*W  "J  *.y*J  "  P*rmW*r*3m**T*jrV*m*«yrmymy*y  If  g»r*  V*  >*  r»V»  y»  vy  If  tf  yjf Jf .tf Jf  tr»yi 


\ 


2)  In  the  same  v ay  which  was  the  case  with  the  forming  of  a  model 
for  steady  state  combustion,  in  the  combustion  of  propellants,  the 
phenomenon  of  molten  binder  agent  covering  oxidizer  surface  area  should 
also  be  considered  in  the  forming  of  a  model  for  unstable  combustion  of 
solid  propellants. 

3)  Using  the  unusually  regular  combustion  of  the  congealed  phase 
reaction  and  opposed  gasification  covered  by  molten  binder  agent  to 
replace  localized  shut  down  is  a  precise  and  necessary  method.  This  is 
true  not  only  for  steady  states,  but  for  non-steady  states  as  well. 

4)  The  model  presented  in  this  article  is  still  awaiting  further 
experimental  testing  and  perfecting. 
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SUMMARY 


In  this  article,  under  simulated  engine  combustion  chamber  p  ',e--3  -3  ar*  -5 
conditions,  we  carry  out  an  analysis  of  and  calculations  for  the  dynuiios 
of  the  omni-axial  rotation  of  flexible  jet  nozzles.  The  important 
contents  include:  analytical  calculations  of  measured  data;  analytical 
nulculations  for  the  center  of  instantaneous  oscillation  of  jet  nozzles; 
and,  calculations  of  the  moments  of  force  of  oscillation.  Besides  this, 
we  also  calculated  the  paraaeters  that  follow:  the  angle  of  oscillation 
and  angular  velocity  of  jet  nozzles;  azimuth  angles  of  oscillation;  tie 
length  of  operating  tubes;  eccentricity  of  thrust;  and,  axial  and  radial 
displacement  of  jet  nozzles  under  differing  pressures  of  contents,  etc. 

This  method  of  calculation,  after  having  appropriately  resolved  the 
p ” sbl e a  of  installing  telemetry  pickups,  is  still  suitable  for  use  with 
hot  test  beds  and  is  also  suitable  for  calculations  and  measurements 
under  conditions  of  static  oscillation  of  flexible  jet  nozzles. 

1,  Calculations  of  Adjustments  in  the  Amounts  of  f)  isplace  a ont  and  Angles 
of  Oscillation 


(1)  Selection  of  Measuring  Equipment  and  Coordinate  System 
The  apparatus  for  measuring  simulated  axial  oscillation  of  flexible 
jot  nozzles  is  -is  shown  (?ig.  1  ) .  The  flexible  connector  head  connects 

"  e -sp  jo  ti  vely  with  the  end  ring  and  the  lower  flat  plate.  They  and  the 

container  shell  body  together  fora  a  sealed  high  pressure  vessel.  To  the 
hot  to  a  flat  plate  is  firmly  attached  a  rod  (simulating  the  jet  nozzl 0 . 

Its  axis  line  and  the  axis  line  of  the  flexible  connector  head  are 

congruent.  Its  point  of  contact  on  the  Z  plate  is  the  origin  of  the 


content  pressure 


zer o . 


The  axis  OZ  is  then  congruent  with  the  axis  line  of  the  rod. 

The  top  end  of  the  rod  and  the  axis  line  are  installed  at 
•’ '  gh  f  tngl  •?  •:  f  i  th  *  "1  a  t  plate,  «’h  i  j  n  '  ?  called  the  Z  plate.  The  axial 
1  ■  n  ■>  )f  th  t  lisplacenent  sensor  z,  is  congruent  witn  the  axis  line  of 
th  e  od  .  "”1  =>  direct-"  on  end  amount  of  thrust  are  cont^oll-d  by  actuator 

tubes  '  an  1  2,  which  are  r  > spec ti vely  positioned  in  plane  3  ox  z 

ani  ny 7  •  I n  this  way  a  fixed  coordinate  systei  O. Y  YZ 


is  then  established.  Two  flat  plates  are  installed  on  the  rod. 

They  are  respectively  perpendicular  to  the  X  axis  (called  the  X  plate; 
and  the  Y  axis  (called  the  Y  plate).  Displace  neat  sensors 
and  Yj.K,  are  respectively  perpendicular  to  the  X 

plat-  and  the  ^  plate.  Moreover,  lines  extending  trot  A, 

and  ,  and  A',  and  V,  connect  respectively  with  two  points  on 

the  7,  axis.  Their  points  of  i ntersect ion ,  with  coordinates  in  these 
dimensions,  are  (0,0,-Z,),(0,0,-Z,)  .  An  extension  of  the 

axis  line  of  the  displacement  sensor  intersects  with  the  0^} 
coordinate  plane  at  point  (O,  —  Y,,  —  Z,)  .  The  distance  fro. a  the 

plate  and  the  Y  plate  to  the  coordinate  point  is  F.  in  both 
cases.  The  other  two  displacement  sensors  H>-  y.nd  respectively 

measure  the  two  actuator  tubes  fixed  to  rod  points  Q,  arm  0:  in 

terns  of  the  ancmt  of  axial  displacement  caused  by  elastic  deformations 
of  t ne  vessel . 

Besides  this,  the  radius  of  measurement  contact  heads  for 
displacement  sensors  X ,*X ,  is  in  all  cases.  The 

radius  of  the  measurement  contact  head  for  sensor  Z,  may  be  selected  a 
will.  Let  this  value  be  ft,  .  The  rest  of  the  ins  -.alia Mon  constants 
and  structural  constants  are  a;  shown  in  the  Figure. 

F,  (or  p,  (,  and  F,  (or  jD>  )  are  activating  forces  (or 
pressure  differential)  sensors. 

(?)  Analysis  of  D is pi ace  nent  Sensor  Measurement  Data 

Amon-r  th J  eight  displacement  sensors,  except  for  the  measurement 
•1  tta  fr  i  !  H,  n ’  M,  r.ii oh  there  is  no  need  to  correct,  the  measured 

1-iti  >•>;.'  ml  1  >o(i  the  other  six  sensors  certainly  does  not  perfectly 
reflect  the  aaount  of  displacement  in  the  oscillation  of  the  rod.  This 
is  due  to  the  reasons  set  out  below. 

(1)  The  measuring  contact  head  of  the  sensors  is  a  hemisphere. 

Ac  lolling  t )  changes  in  the  angle  of  oscillation,  contact  points 
experi-Jic1  displacement  along  the  hemisphere,  introducing  an  additional 
amount  of  displacement. 


Fig.  1  Schematic  Diagram  of  Test  Equipment  1.  Plate  2.  Rod  3.  Lower 
Flat  Plate  4»  Upper  Compression  Ring  5*  Flexible  Contact  Head  6.  High 
Pressure  Container  Shell  7.  Actuator  Tube  2  for  Controlling  the 
Direction  and  Amount  of  Thrust  8.  Plate  9.  Diagram  10.  Plate  11 . 
Actuator  Tube  1  for  Controlling  the  Direction  and  Amount  of  Thrust 


(2)  The  measurement  contact  heads  for  the  various  displacement 

sensors  Xf>Xi%X,xYf>Y,  and  the  contact  points  of  the 

experimental  mea.«nr  ?  .  ent  plates  are  at  a  distance  from  contact  points  on 
the  rod  axis  line  X ,  and  from  the  symmetrical  surface  of  the  rod  so  that 
the  distance  varies  with  changes  in  the  angle  of  oscillation. 

(3)  Because  of  the  fact  that  the  axis  lines  of  X,  Y,  and  Z  rotate 
at  the  same  time,  their  mutual  influences  cause  plates  X  and  Y  to 
correspond  to  the  slant  which  occurs  in  the  sensors.  This  introduces  an 
amount  of  non-displace aent  in  th  i  s  ensors. 


i 

i 
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(3)  Adjustment  Calculations  for  Angle  of  Oscillation  and  Amount  of 
Displacement 

Let  the  angles  of  rotation  for  the  axis  lines  X,  Y,  and  Z  be  0*a*$ 

.  In  order  to  simplify  analysis  and  calculations,  we  <iH<e  use  of 
the  geometrical  principle  of  superposition  to  make  a  two  step 
calculation.  After  that,  we  superimpose.  That  is,  first  we  make 
calculations  for  the  situation  in  which  a*0  and  0*0  ,  withfl«*0 

Then,  we  make  the  calculations  for  the  situation  in  which 
a* 0,0*0  and  8*0 

(1)  Adjustment  Calculations  for  the  Situation  in  Which  c#o,/9#o 
and  we  assume  0=0 

There  are  four  types  of  factors  in  omni-axial  oscillation  which 
cause  adjustments  to  the  displacements  measured  by  the  displacement 
sensors:  (a)c>0,/J<0»(b)a<0,  /?<Oj(c)a<O,0>O|(d)a>O,/J>O. 

To '’■aula  derivation  brings  up  one  type  of  situation. 

Eccentric  rotation  of  the  axis  line  within  the  first  apparent  limit 
a>  O,0<O. 

If  we  assume  that  the  cross  section  of  a  rod  is  a  square  23.  on  a 
side,  then,  Fig.  2  is  a  cubic  diagram  of  angle  a  t  id  angle  0  ef  the 
deflection.  Fig-  3  is  a  horizontal  projection  diagram  corresponding  to 


Definition:  The  plane  which  holds  the  displacement  sensors  Xt 

and  Y,  (or  X ,  and  Y,  )  is  the  horizontal  plane  of  test 
measurements.  The  same  rod  cross  section  used  in  it  is  reprsented  by  uEFGH 
If  we  assu'ie  that  axis  line  Y  rotates  through  angle  a  and 
axis  line  X  then  rotates  through  angle  P  ,  the  change  in  the  horizontal  (^1 
plane  of  test  measurements  is:  squar°  DE,F,G,Ht~  elongated  square 
zsE'F'G'H quadrilateral  &EFGH  .  This  corresponds  to  a 
change  in  the  contact  point  of  the  experimental  measurement  head  of  the  X 
displacement  sensor:  At-*At-*Al0 

At  point  a,  •  a  =  0,0  =  0  •  X  sensor  measurements 

give  us  a  displacement  amount  A'.  =  0 ,A,0  =  R 

At  point  A,.  a> 0,0  =  0  ,  The  neccrb  amount 

disolacement  is  X'.  The  actual  amount  of  d  i  splae  e  lent  Is  A' '  =  0,0 


The  men'ur-^d  amount  oJ- 


One  of  the  additional  amounts  of  di splace  aenl 


A,0,~  AtO= R/oosa-  R 


Pro  (i  the  diagrams  it  is  possible  to  see  four  amounts  of  additional 
displacement  caused  by  the  displacement  of  the  contact  points  on  the 
contact  heads:  A,At  =0'A;  -0'A,  =  ft,/cosa  —R, 

X'  =  X  +  (ft+ ft,)(l/ccsa  —  1)  (1) 

At  point  Ai.  :  a>O,0< o  .  Let  the  numerical  data  measured 

by  the  sensor  be  X"  .  -From  Fig.  2  or  Fig.  3  it  is  easy  to  see  that 
the  value  of  the  difference  between  X*  and  X'  is  A'iAt, 

These  are  three  of  the  amounts  of  additional  displacement  caused  by  angle  /?. 
on  sensor  X  .  From  Fig.  2  we  obtain  the  set  of  relationships 
set  out  below; 

From  RtAA,A',A',  we  get:  '  A',A‘,‘  —  A,A',tg0  =  Ytg/3 

From  fttAA,A’,At  we.  get:  A'lA,  =  A'tA'itga  =  Yigatg/3  (p) 

Fro  a  RtAA,A,A'i  we  get:  AtA,=A,A,tgr  =  Ytgr  (3) 

Compared  with  equation  (2),  equation  (3)  has:  r  =  arctg [tga  •  tg0]  (4) 

Because  A\A,=Yigr<0  >  We  also  give  consideration  to 
satisfying  the  geometrically  equivalent  relationship,  and  we  then  have: 

'  X*-A':At~X'  (5) 

After  we  now  do  some  more  analysis  of  the  deflection  angle  , 

the  four  additional  amounts  of  displacement  which  are  induced  in  the 
displacements  of  the  contact  points  on  the  measuring  contact  head  of 
displacement  sensor  X  can  be  seen  from  Fig.  4.  The  contact  point  moves 
f  -on  point  a,.  to  point  A ,,  :  that  is,  it  rotates  through  angle  X.A'-0,lA,  =  r,0KAl 
=  ft, cos  a,  Because  of  this,  there  is  always  a  positive 
mount  0"  increase: 

J?,oosa(tyeosr,l?-l)  '  (6) 

When  one  combines  this  with  equations  (l ) (2) (5) (4) (6)  one  gets: 

Below,  when  we  make  adjustment  calculations  for  the  amount  of 
dlspl  '..cement  measured  by  the  Y  direction  displacement  sensor.  We 
similarly  have: 


At  point 


F= 0  5.0  =  /? 

Y'-Y+iR+Rt)(  1 

'  cos  ft  ) 

At  point  Bt  :  BjB,^Xtgr 

■  Xtg  r<o  ,  it  satisfies  the  equivalent  relationship 

Y"  +  575,=  Y' 

Y'-Xig  r  =  Y' 

The  aid i t ional  amount  of  displacement  caused  by  the  angle  a  is: 

R,cos0(  1/cosa—  1)  (10) 

From  the  various  equations  (3) (9) (10)  one  then  obtains: 

Y inY*  —  X  igr  —  (R+ R,)(  l/cos/3—  1)  —  R, cos  0(1 /cos  r  —  1 )  0  0 

(2)  Calculations  to  Adjust  the  Amount  of  Displacement  When  the 
Angles  a*0*O  A re  Simultaneously  Wot  Zero 

As  far  as  an  ideally  flexible  contact  head  is  conerned,  one  has 
only  the  times  when  the  azimuth  angle  of  oscillation  0.  is  45  degrees 
and  225  degrees  ,  making  use  of  forces,  the  combined  force  of  which 
crosses  over  the  axis  of  the  jet  tube.  In  the  case  of  oscillations  with 
other  azimuth  angles,  making  use  of  forces  the  combined  force  of  which 

not  at  all  cross  the  axis  of  the  jet  tube , adds  an  additional  moment 
of  force  to  the  axis  line  of  the  jet  tube.  In  this  way, the  jet  tube 
rotates  around  the  Z  axis  line,  the  angle  of  which  is  designated  0 
During  omni-axial  oscillation,  it  is  only  possible  for  the  angle  8  to 
be  positive  or  negative.  If  one  gives  consideration  to  a.0  ,  the 

overall  influencing  factors  can  be  grouped  into  six  types  of  situations, 
that  is, 

0“  =  tf,<45“  6<O,(a>O,0<O) 

45°<it,<90e  0>O,(a>O,5<O) 

90°<yi,«Sl80o  0>O,(a<O,0<O) 

180°  <  tp,  ^225°  6>  0,(a<0,5>0) 

225°  <^*5270°  0<O,(a<O,5>O) 

270“  <  \0,^36O“  6<O,(a>O,0>O) 


(8) 
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(9) 


Pig.  4  Movement  of  Contact  Point  on  Measuring  Contact  Head 

Pig.  5  Diagram  of  Rotation  of  Horizontal  Test  Measuring  Surface  «<o 

Doing  through  a  practical  derivation  clearly  shows  that  the  six 
situations  all  possess  the  same  type  of  expression.  In  this  case, it  is 
only  necessary  to  do  calculations  for  one  situation. 

0°^^,<45°,a>0,/3<0,e<0,A'>0,K>0 

0  is  the  angle  of  rotation  around  the  OZ  coordinate  axis.  It  is 
o  r.  Vi  it  th  ho  ’ i  cental  test  aeasurement  surface  rotates  through 
ingle  0  in  the  same  way.  Analysis  of  the  front  surface  down  to  the 
mutual  influences  of  the  angles  of  deflection  “  and  0  cause  the 
horizontal  test  measurement  surface  to  give  rise  to  deform  it  ion.  After 
maxing  adjustments,  a  and  0  are  used  independently  to  cause  the 
horizontal  test  measurement  surface  to  give  rise  to  deformations  in  Oj E'F'G'H" 
as  in  the  diagram,  and  rotate  it  through  the  angle  0  ; 

f’n.it  i  s ,  becoming  zjE"'F'"G'"H'"  ,  as  shown  in  Pig.  5- 

proa  ^ig.  5  one  can  see  that  the  angle  6  causes  the  X  sensor  and 
v  -?  e  'i  s  >r  t  •>  r  ^  3  o  aot  i  vel  acnnire  an  additional  amount  of  displacement 

4/  U  X 

been  ling  A'iA,  and  jJ\B 
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Solving  for  A'iAu  t  the  D  point  coordinates  are 

(X0—X  +OtDcoeO 
W  d  =  Y  —  0,Dsin6 

A.O.-O.D,,  X.4;0  =  Jf'  =  A’+(7?  +  /?.)(— l - lW'+^.O, 

\  cos  a  ) 


:.0,D=A,0,  =  (R+Rt)(  — ]-  -  -l) 

\  cos  a  / 
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,  also  because 
therefore,  the  D  point  coordinate  should  be: 


DC«0,D  sin0<O  , 


'  XD=X+  +  -i\cos<> 

I  \  cos  a  / 

* 


(12) 


Make 


and  H'E "  cross  at  point  J.  Then,  the  equation 


for  the  straight  line  DJ  is: 


K__ii 


hlso  one  anr 


X,=  X'  =  X  +  (R+R,)(— i - 1) 

\  cos  a  / 


(13) 


(H) 


Take  equation  (14)  ani  substitute  in  equation  (13): 

*— id  (/?  +  /?,)(— ^—-l)(l-cos0V  +Y  +  (R  +  R,)(  - 1 ) sin i 

tgPL  \  cos  a  /  j  \  cos  a  / 


(15) 


wro'H  Rt.\A"A,J  we  get: 


A;A,  =  YjigO 


(16) 


,  "H 
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Prom  Pig.  5  we  get  the  geometrical  relationship 


a,o+a‘;a3=a;o=x' 

A,0=X%A\'A,  =  y,  tg$<0  :.A,0-A:A,-A';A:  =  A,0 


Substituting  the  various  algebraic  quantities  and  simplifying,  we  get: 


*  =  «gr+,g0)'-(tf+/o(  1  -l  )!-/?,«*«(  JL_j\ 

\  cosa  /cos#  leas  r  / 


(17) 


In  solving  "  >r  ,  we  use  a  process  identical  to  that  use!  in 

solving  for  A- A,  ,  and  we  obtain  an  S  point  coordinate  which  is: 


(18) 


Let  E'"F'"  and  E0F"  cross  at  point  Q.  The  straight  line  equation 
of  3Q  is 

Y~tgO(X-Xs)+Y,  (19) 

The  0  point  coordinate  is:  (20) 

T  it  ?  ^nations  (20)  (id)  and  substitute  in  equation  (19): 


)  a 


Rt&B,B'lQ 


is  possible  to  obtain  B'iB,=  Xoig6 


(22) 


The  geometrical  relationship  has  the  form:  B,0- B\B,=  B'-O 
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s  V 


•.■B"B,-XoigO<0 


;  at  the  same  time,  after  giving  consideration  to 
the  superposi tion  quantities  for  B'iB,  ,  one  then  has: 


(23) 


Equations  (17)  and  (23)  are  the  equations  for  calculating  the 
amount  of  displacement  adjustment  which  is  measured  by  the  displacement 
sensors  in  the  X  and  Y  directions  during  omni-axial  oscillation.  The  two 
are  both  functions,  and,  set  up  together,  they  can  be  solved  for  an 
equation  through  which  the  amount  of  displacement  (X ,K)  can  be 
calculated  from  the  measured  amount  : 


X  n  = 


X’r,-Y)',{lgrt+ig6t)-  (R  +  Rt)'\  (-JL - _!) 

<x*s9,i.  \  oosa,  )  \  cos/9,  / 


1-ig*;* 

(tgr,  +  tg0,)J  i  )[C06a,- cos/9, (tgr,+  tg0.)] 


(24) 


+  ig’9i 


Y’,:-X::(tgrl-igO.)-(R  +  Rt)  *  r  (  *  -lW-  -i - 1 ) 

y  CO s9.  \  cos/3,  '  \  cosa,  ' 


1-tgV, 


(25) 


(tgr.  -  tg0,)] - l)[oos0,-cosa.(igrt- tg0,)] 


+  tg  '6i 

.a  the  equations  /=  l»2,i  —  0,1- 


a.,/3,,r,,0, 


are  calculated  from  the 


RiR* -  constants; 

formulas  given  below. 

(3)  Calculation  of  Angle  of  Oscillation 

Prom  the  definitions  of  the  various  angles  of  rotation  one  can 

have : 

laa  =x-‘lx  x, .  =  x;:-x;;-(K;;-K;:)(tga,tg^.  +  tg0,) 

L  Zx—Zx  (1  —  tg’c.tg^.  +  tg'flJZ,  (26) 


a  =  _  Xu-_Y_l  =  Y'x'-Yy,+(X'x:-X',’,)(iga,ige,-ig9,) 
Zx—Zx  U  —  tg'a.tg'/J.  +  tg^.OZc 


(27) 


Calculation  of  the  Angle 


o, 


Prom  equation  (17)  it  is 


possible  to  obtain: 


A’.cosajf  — - l) 

\  COST,  ' 


-^2i*A';:-Kli(tgr1  +  tg0i)-(/?  +  /?.)( — - - 1)— 4— 

V  COSOi  )  COS0i 


COS0, 


Through  deduction  it  is  possible  to  obtain: 


x  u=x  \ : —(^h+k,  )(tgr,+tg  0, )—'/?+/?,)( — - — iV— - 

\  cos  a,  /  cos  0, 


—  /?,cosa(^ 


1 


oosr  j 


(28) 


?coi  the  i e f L a  1 1 i. o a : 


ig^i 


xv.-xw  l 


2K. 


-2l8ff 


(29) 


Equations  (26) (27) (29)  are  a  set  of  implicit  transcendental 
equations.  Using  the  iterative  substitution  method  in  a  computer 
application,  one  solves  for  values  of  .  After  this,  one 

then  substitutes  equations  (24)  and  (25)*  It  is  then  possible  to  solve 
for  X,i  and  K,<  for  the  various  amounts  of  displacement. 

(4)  Adjustment  of  the  Numerical  Data  Measured  By  the  Displacement 
Sensor  Z, 

Definition:  The  angle  of  oscillation  is  the  angle  included  between 

the  geometrical  center  line  of  the  jet  tube  and  the  0 Z  coordinate  axis. 

It  is  expressed  by  the  use  of  9 

The  numerical  data  recorded  by  the  Z,  sensor,  alone,  due  to  the 
-ngl  -?  of  slant  on  the  Z  plate,  introduces  an  additional  aQount  of 
1.  i  spin  ce  dent .  Moreover,  this  angle  of  slant  is  always  equal  to 
oscillation  -njle  $  -  In  order  to  solve  for  this  oscillation  angle, 

o  i e  aus  t  fi’st  solve  for  the  direction  number  of  the  geometrical  center 


Oscillation  angle  t 


COSO,  = 


*JA\  +B]+C' 


(31  ) 


Adjustment  calculation  for  Z\\  ;  Zu-  Z\\—  k{ — - - (32) 

\  coso,  ) 

In  a  heat  test  bed  situation,  it  is  possible  to  make  use  of  two 
dfsplaoeaent  sensors  Z%*  -uni  Z,»  la  order  to  replace  sensor  Z»  .  < 
ani  Z,,  are  installed  symmetrically  on  the  two  sides  of  the  jet  tube  in 
order  to  avoid  ignition  flame  heat  corrosion. 

As  "Pig.  6  shows 


'  T-S 
a  ni 


->  v  a 


' " '  Chart  of  Sensor  Distribution  (b)  A  Comparison  of  Sensors 
to  Sensor  Z, 


It  is  easy  to  demonstrate  that,  between 


■I  4 


z 


and  Z, 


there  is  the  relationship  shown  below: 

z;'=(z;;+z;;)/2 

(5)  Adjustment  of  the  Coordinates  of  Actuator  Tube  Piling 
Cjnt'ol  the  Direction  and  Amount  of  Thrust 

Due  to  the  ’esults  of  pressure  contents  and  the  forces  us 
actuator  tubes,  the  displacement  of  the  rod  filing  po:nt  along 
r  e  a cues  ab  o v s  5  m l 


(33) 

Points 


This  not  only  produces  an  influence  on  the 


p  r  o  0  e 


the  actuator  t  foes .  It  is  also  related  to  the  ii-ec 


forces  used  and  the  force  momentum.  Therefore,  it 

necessary  to  make  adjustments  to  the  measured  amounts.  Host  of 


ed  with 
the  a  a  is 
shock 
ion  of  t 

is 

the 


amount  of  the  change  is  along  the  11  faction  of  the  axis.  Therefore,  the 
adjustment  is  approximately: 


(ZQii  —  Hq  +  H\i 
'■Zq,  (  =  HQ  +  H’ii 


(34) 


In  the  equations,  Hq  is  the  coordinate  of  the  rod  fixing  point. 
are  the  values  of  the  amounts  measured. 

(6)  Calculations  of  the  Motive  Forces  and  F ,<  TJsad  to 

Control  the  Amount  and  Direction  o f  Thrust  F, 

As  motive  force  one  can  choose  to  use  the  tensile  pressure  sensor 
contact  measurements  obtained  as  F  ,t  and  F ,i  .  Or,  one  can  measure  the 
pressure  d i f f er ent ial  P  on  the  two  sides  of  the  piston.  The  force  used 


is  then: 


Fn  =  SPi 


In  the  equation,  S  is  the  area  of  the  piston,  ;  =  l,2,i=  1,2,-  n. 

TJp  to  the  present  time,  we  have  solved  for  the  following  parameters. 

Zq  1 1  \Zoti  %}*  it\F  n\0 1  \  ft  t  1  \  f  t  \  i  t  %At  %  lit 

and  C  ,  as  well  as  constants  to  be  used  in  later  calculations.  It  is 
not  necessary  to  explain  these  for  our  purposes  here. 

2,  Calculation  of  Center  of  Oscillation  for  Flexible  Jet  Tube 


Tl-C  OTiOBPf  OP  A  CCITCI  )  P  030  [hi.  AT  [04 

As  far  as  the  concept  of  a  center  of  oscillation  for  flexible  jet 
tubes  is  concerned,  both  inside  and  outside  China,  there  have  been 
treat  i  mbs  in  a  amber  of  technical  reports.  But,  these  reports  have  not 
found  a  clear  and  precise  definition.  Moreover,  there  have  been  several 
methods  put  forward  which  are  obscure  and  unclear.  Before  making 
c  ileal at  ions  of  the  center  of  oscillation,  we  plan  to  have  a  clear  and 


precise  concept  to  act  as  a  theorm. 
calculations . 


a  i  -ta  .  i- 


‘or  analysis  and 


-M  t  J  l'\  1'  Ci 


ad er stood  in  the 


following  terms.  In  the  case  e *  r i g : 1  lediss,  there  are,  of  course  plane 
oscillations  taking  place  as  well  a'  \  it',  al  is-'illations  (omni-axial 
:>  ec  illation)  •  Both  of  th-se  exist  simply  as  a  ser :  es  of  instantaneous 
a x i. al  rotations,  and  no  i.  os  tant  snmus  center  of  oscillation  exists. 
Obviously,  a  jet  tube  is  a  rigid  body,  and  its  center  of  oscillation  can 
ciianqe . 
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Flexible  jet  tubes  also  have  a  point  which  functions  as  their 
theoretical  center  of  oscillation.  If  the  design  and  nanufacture  of  the 
flexible  connector  heads  are  both  absolutely  ideal,  this  leads  to  the 
flexible  jet  tube  rotating  arounl  the  same  fixed  point.  What  are  called 
its  instantaneous  axes  of  rotation  all  pass  through  its  theoretical 
center  of  oscillation. 

Absolute  idealization  of  flexible  jet  tubes  does  not  exist.  In 
reality,  due  to  content  pressure  (or  combustion  chamber  pressure) 
changes,  flexible  jet  tubes  are  made  to  produce  relatively  large  amounts 
of  axial  instability.  Besides  this,  flexible  connector  heads  do  not  have 
uniform  distribution  of  material  around  their  circumference,  and,  during 
oscillation,  their  acceptance  of  forces  will  not  be  entirely  symmetrical. 

It  is  also  impossible  that  the  deformations  in  the  flexible  connector 
heads,  which  a- e  caused  by  this,  should  be  symmetrical  and  uniform.  This 
will  also  cause  the  instantaneous  axis  of  rotation  of  flexible  jet  tubes 
not  to  pass  through  their  theoretical  center  of  oscillation.  When 
flexible  jet  tubes  undergo  omni-ax  Lai  oscillation  under  various  types  of 
content  pressures,  the  series  of  instantaneous  axes  of  rotation  will 
form  a  cylinder-shaped  included  figure.  This  figure  will  be  symmetrical 
around  the  center  of  oscillation  which  will  be  included  in  it.  This  is 
called  the  "oscillation  center"  for  short.  This  is  also  nothing  else 
except  the  actual  oscillation  center  of  the  jet  tube,  the  size  of  the 
geomet r Leal  dimensions  of  which  reflect  the  qual  ity  characteristics  and 
rigidity  of  the  flexible  connector  heads. 

If  one  is  considering  directly  solving  for  the  geometrical  range  of 
the  included  center  of  oscillation,  it  is  necessary  to  have  dense 
instantaneous  axes  of  rotation.  It  is  also  necessary  to  have  large 
amounts  of  experimental  data.  It  is  difficult  to  get  to  this  point,  and 
it  is  also  not  economical.  On  the  basis  of  this  special  characteristic  of 
jet  t  ib  is  -otating  around  these  fixed  points,  it  is  possible  to  reach  the 
conclusion  which  follows.  The  i nstantaneous  oscillation  center  of  the 
ge  d  ae  t  r ;  cal  center  line  of  jet  tubes,  is  always  a  reduplication  of  this 
fixe!  point.  This  ^  i  x-d  point  what  is  called  the  oscillation  center 


of  this,  we  then  use 


i nstantaneous 


enter 


geometrical  center  line  of 


jet  tube  to 


replace  the  point  center  of  oscillation  of  the  jet  tube.  By  means  of 


this  method,  we  use  a  limited  amount  of  experimental  data  to  solve  for 
the  geometrical  size  range  of  the  jet  tube  oscillation  center.  Moreover, 
this  supplies  a  reference  point  for  solving  for  the  L nstantaneous  moment 
of  force  of  oscillation. 

During  oscillation,  if  the  geometrical  center  l lues  cross  each 
other,  the  results  of  the  calculations  discussed  above  are  co aoletely 
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accurate.  However,  it  is  also  possible  for  them  to  be  mutually  off  and 
not  crocs.  In  such  a  case,  one  then  taxes  their  "average  oscillation 
plane"  and  solves  for  the  oscillation  center.  Obviously,  the  results  of 
this  are  approximate.  The  degree  of  this  approximation  varies  with  and 
increases  as  the  degree  of  graduation  of  the  measured  angle  o  f 
oscillation  is  reduced.  Experimentation  el  ef’ly  Vw  ■  that,  vnen  the 
measured  angle  of  oscillation  is  changed  f -'on  ‘5  degrees  to  1  .5  degrees, 
the  distance  between  two  adjacent  geomet r i  ;  tl  :>oter  lines  drops  from 
0.5 -nm  to  0 . 2mm , roughly  speahing. 

CALCULATION  OP  THE  OSCILLATION  CENTER  OP  PuEiiti.E  -IE1  TUBES 

In  order  to  simplify  the  discussion,  let  to-  geo  net r leal  ;  enter 
line  of  the  jet  tube  be  the  Z'  axis.  This  is  riot  the  saae  as  the 
axis  of  the  coordinate  system. 


The 


axis  equation  is: 


X-Xu  _Y-Y,±Z-Z 


A, 


B, 


(35) 


’he  equation  for  the  Z  plate  plane  is  A,X  +  B,Y  +  C(Z  —  Z,,)  =  0 


The  set  of  simultaneous  equations  (35)  (3b)  gives  us 
coordinates  for  the  0,  point.  The  0.  p  ’  ’  i  ’  ‘  n  -.t  ? 


,  _CZ,i-AtQi-BiPl 


A,  K  i  +  p 1 1 1  (7 

x  u ,  ~  A  ,  r , ,  Q , 


O’  • ,  —  L  ^  ( +  A* , 
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kv. 


5o  ! 


Vc'v  '  /O-  V-  -  -  - .  ■ 


.‘jV1 


iv--', 
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In  the  equations 


K,=At/CiQl  =  X>,-K,iZ^Li  =  Bi/C%Pi  =  Yu-l.Z> 


The  0<  point  is  a  point  on  the  Z'  axis.  In  the  initial 
conf iguration,  it  is  a  duplication  of  the  origin  point  of  the  coordinate 


s ystem. 


If  we  postulate  two  mutually  adjacent  instants,  the  plane  at  0  i  which 


contains  the  instantaneous  center  of  oscillation  is  n„ 


We  then 


nave : 


(jc — x,<)’  +  (.y—.y..  )’  +  (*— *.<)'  =  (*— 

After  r ear rangement ,  we  gst:  *><«  A,iX+B,iy+C,iz= D,t 


In  the  equations 


(Ail  —  Xfl  Xu— if  Bu  —  yoi  y*i—  l  1  Cli  ~  *»i  Z*i  - 
(D,.  =  (jc, )+y»)  +  z,\  -xt]-,-y,)-i-z.)-,)/ 2 


In  order  to  solve  for  the  i ns tantaneous  axis  of  rotation  of  jet 
tubes,  one  must  still  find  a  second  point  on  the  plane  on  which  we 
find  the  center  of  oscillation  for  the  time  interval  during  the  same 
period.  Let  this  be  the  intersection  line  of  the  plane  *»/,*»<  and  the 
plane  *.<  ,  and  one  then  has  the  instantaneous  axis  of  rotation. 

Because  of  the  fact  that  a  measurement  error  in  the  system  is 
i.inavoilabl  ^  it  is  only  when  n,i±7i„  that  the  degree  of  precision 

associated  with  the  axis  of  rotation  which  we  find  is  optimum.  The  O, 
point  is  a  Pixel  point  on  the  jet  tube  axis  line.  Points  which  satisfy 
the  joilition  will  vary  with  changes  in  the  oscillation  of 

the  jet  tube.  The  method  for  the  fixing  of  coefficients  is  none  other 
than  an  unceasing  pursuit  of  these  points. 

If  we  take  the  plane  ,  and  set  *«  f  nlt 


A,iX  +  B,iy  +  C,,z  =  D, 


If  we  solve  the  simultaneous  equations  (35)  all  [}')) 
j  ? :  ‘  i  *  i  it  >  3  ;tion  point  of  *«•  and  the  Z;_.  a< 
When  we  solve,  we  get  the  result: 


E„(xt,,y„,z„ 


’»  _**,*«  .*»  ,*• 


The  coordinates  for  point  £»•  : 


in  the  equations 


'  _  _  C Dti  —  A i —  B i , P , )  _  is  n  r\ 

2“  AuK +  B,,Li-,  +  C,i  K"  *'  0,1 

x, ,  =  K  ,,D,  +  Q,, 

y, i  —  K.  u  D.  +  Qu 


K  t  <  +  C,<) 

K  K  h  Ku  =  Li-,Ktl 

On  *■  K i-,Q„  -r  O.-i  Q».  =  L.-.Qu  +  /J, 


M ■  ate e  the  *  jtion  point  of  plane  and  axis  Z’,  be  £».(*».  »>».**#.)# 

Whin  n  nl/j,  y;  ^et  th  -*  "^sult: 


Coordinates  for  point  E ,, 


f:  _  d..-a,.q,-b,.p. 

”  A,,K ,  ~  B UL,  +  C,, 
h,.*L...D..  +  P., 


in  the  equations: 


A,,K,  ■+■  5 1 ( L t  -t-  C , , 

1 #  “  ^  (liti 

f-  ii  e  i- 1  i- 1 1  i 


o.o-r* 


p  -  _-Ld'.o, j 
1  '•  'A„K,*H..L,*C„ 

P„-Q,  +  K.P„ 
P„~P,  +  L,P,. 


K..L..P..U. 


• .  *  •  ■»•»  y  , 
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4fte-  o'H  takes  o  quations*  (40)  and  (41)  and  substitutes  them  in  an 
expansion  arrangement,  one  has 


in  the  enuatlons 


UtlDm\  +  VttDmt  +  Wti  =  0 

U„ -/,*+/,*+/,!  -A',; -K,]- K,\ 

V h  ~  2[/Ci ,( x(l_ ,  0 , , )  +  Q.,)  +  /e>,(*„-, —  On) 

■£•!<(•*#<  Pn')  LttCyn  ^ij) 

C y$i  -  i”~Ojt),“~(^i (—  i ^ * Oi i ) * 


)l  /i.  ng  equation  (42)  we  have: 


2U,t 


(43) 


O,,,..  ig  then  the  fixed  coefficient  we  were  solving  for.  There 

;v)  wts  of  solutions.  However,  cor  responding  to  our  actual  problem, 

to-*  --elution  -ef  in  which  r„  and  *».  are  simultaneously  smaller 
sere  is  aooi  i  ruble  Therefore,  we  take  the  two  solution  sets  of 
,  vr  -  ~  i  a  a  1 1  3  lb  s  t  L  tuts  then  respectively  into  equations 

uni  4‘).  ?i‘oa  tue  principles  discussed  above,  we  decide  and  make  a 

■  ‘-.on.  Because  of  this,  after  we  make  a  decision,  the  solutions  of 
i-:1  r  •;  h ’  :h  -i  •  ;h  )  ;  -> ,  then  auks  it  possible  to  solve  for  the 
n,.,  on  which  is  located  the  ;s  i  iul taneous  oscillation  center  for 


Moreover , 


*»,  JLw,., 


(  x  -*..)’  +  (  y  -y„y  a  (t  -t„  V  -  (  * -x„y  +  (  y  -  _v„  Y  +  (  *'  -  )’ 


At,x  +  Bt,y^Ct,y  Dt, 


SOT 


in  the  equations 


(Atl  =  xu  —  ,  B„  *=  y%,  —  ytl ,  C,(  =*  ztl  —  z,, 

iD„  =  (x,l  +  y,l  +2,)  —  jc.  f  —  y,;—  z,])/2 

The' line  of  intersection  between  planes  -*u  and  is  the 

instantaneous  axis  of  rotation  we  are  solving  for.  Below,  we  also  solve 
for  the  plane  of  oscillation  of  the  geometrical  center  line  Z'  axis  of 
the  jet  tube,  and  we  let  this  plane  be 

During  omni-axial  oscillation  of  flexible  jet  tubes,  the  most 
general  situation  for  the  track  of  the  Z'  axis  is  a  curved  surface  in 
space.  If  we  take  this  curved  surface  and  divide  it  into  »  partial 
curved  surfaces,  and,  respectively  use  a  plane  to  approximately  replace 
their  various  oscillation  planes  an.  ,  it  is  obvious  that,  when  the 
dissected  areas  n  are  increased  in  number,  the  level  of  the 
approximation  is  raised.  We  select  the  number  «  so  as  to  satisfy  the 
requirement  for  precision  which  happens  to  exist.  In  order  to  solve  for 
the  oscillation  plane  n„  ,  one  first  does  the  calculations  set  out 
below: 

(i)  The  direction  number  of  the  common  perpendicular  line 
between  the  jet  tube  axis  line  Z!_,  and  Z\  for  two  adjacent 

i ns  tants . 

Make  the  direction  numbers  of  the  axis  lines  Z\  and  Z'< 
respectively  A rn  and  N u  ,  then: 

N i  =  {A,-,fBt-lfC}  x  {A,,Bl,C)**{A,l,B,i,Cu)  (  45  ) 

in  the  equations 

Au  =  C(B4.,-BlhBu'=C(A/-Af-l)fC,i=Ai-,Bt-AiB,-, 

(ii)  Make  plane  #<->*  contain  axis  ZI-.  .  And,  make  it 

par  allel  to  N,,  Then,  the  equation  for  is 

j  2— X  w-i  t y—Y x— Z, 

Bfm.ll  A  fmm  |,  B  f  mm  I  f  C  =  0 

Am  But  C|( 


•  *  %  -  ’m*  •  „**  * 


«  ,  If  . 


Or  it  is: 


(46) 


A, ,  x  +  B% ,  y  +  C„z  «  D, 


•< 


in  the  equations 


An^Bi-tCti—BuCt  B%t~AuC  Ai-iCti 
C%l-At-,Bu-AuBt-„  Dtt-AtlXu->  +  B>'Y«-'  +  C'‘Z' 


(iii)  Make  the  plane  contain  the  axis  Z '  ,  and  make  it 

parallel  to  N <»  Then,  the  equation  for  R i  is: 
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1  x  X  ,t,  y—Y  lltz— Z , 
R it  ■  A |  ,  B t  r  C 

I  An  ,  Bttt  Ctt 


—  0 


Or 


At,x  +  Bt,  y  +  Ctl2  -  Dti 


(47) 


in  the  equations 


Ah™  B  ,C,,  B,,C  t 
Ctl— AtB  a  AtiBit 


Btt^AuC  AiCu 

Dtl  =  AtlX  ,,  +  B,,Y  H  +  C..Z, 


The  intersection  line  between  the  planes  Ri- 1  and  /?,  are  the 

mutually  perpendicular  lines  ZI-,  and  the  Z\  axis  because  of  the 

fact  that  the  direction  number  of  the  line  of  intersection  is 

equal  to  the  sector  area  (N ,,  x  N „)  of  the  direction  number  of  Z;_ 

and  the  Z\  axis. 

The  various  simultaneous  equations  ( 35 ),(46),( 47 )  necessarily  make  it 
possible  to  solve  for  the  intersection  point  of  the  mutually 
perpendicular  lines  of  the  Z',~,  and  Z\  axes  .  The  reason  for 
this  is  that  the  intersection  line  of  Ri  and  Rt- 1.  is  located  on 
on  the  one  hand,  and,-  therefore,  crosses  the  Z',~ x  axis,  but,  on  the 

other  hand  is  also  located  on  the  plane  R,  ,  and,  therefore,  must 
necessarily  also  cross  Z\  ,  as  is  shown  in  ‘Pig.  7- 

(iv)  In  solving  for  the  intersection  point  1  of  the  mutually 

perpendicular  lines  and  the  Z’,-,  axis,  take  the  equation  for  the  axis 
line  Z;_,  and  change  it  to  be  of  the  form  of  a  parametric  equation. 

Prom  equation  (35)  we  get: 

^x  =  y4,_,/  +  X,,-,  (48) 

j  y  =  £, 

1  *-Cf  +  Z, 
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Pig.  7  The  Mutually  Perpendicular  Lines  of  the  Geometrical  Center  Lines 
for  Two  Adjacent  Instants  Curing  Jet  Tube  Oscillation 

If  we  take  equation  (48)  and  substitute  in  equation  (47),  and,  we  ! 

also  make  t  =  K„  ,  we  then  have: 

K„  =‘(,Dtt—AtiX  i  a- —  i  —  BtjY —  C,tZ,)/(AttAt-i  +  BtiB +  C,,C)  (49) 

If  we  take  and  make  iterative  substitutions  into  equation 

(48)  ,  we  then  obtain  coordinates  for  the  intersection  point  of1  1 

the  mutually  perpendicular  lines:  | 


SXr.i—  i  A,*.tK «i  "t  X  i<  — , 
y  Ei—i**  B  i-,K  t,  +  Y 
2  ft —  C  K ,,  +  Z  i 

(v)  In  solving  for  the  intersection  point  E, 
perpendicular  lines  with  the  Z\  axis: 

Change  the  straight  line  equation  of  the  axis 

form: 


(50) 


of  the  nutually 
Z,  to  a  parametric 


I  x-A,t  +  X„ 
j  y-B,l  +  Y„ 

1  t-Ct  +  Z, 
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If  we  take  equation  (51  )  ami  substitute  in  equation  (46),  and  we 
then  make  t  =  Ktl  f  \ve  have: 


/£■  M  AttX  i  —  CuZ 

"  AtlAt  +  BtlBt  +  CtlC 


(52) 


If  we  take  and  make  iterative  substitutions  into  equation 

(51),  we  get  coordinates  for  the  intersection  point  E ,  of  the  mutually 
perpendicular  lines  and  the  axis  ■ 


=  AtKtl  +  X  u 
yti^  B,Ktt  +  Y  h 

**i =  CKt,  +  2, 


(53) 


It  is  ^asy  to  see  that  the  points  and  E,  are  certainly 

not  the  same  geometrical  points  on  the  axis  line  of  the  jet  tube. 

(vi)  Coordinates  for  the  center  point  E„  of  the  mutually 
perpendicular  lines  E,-.,E, 


Coordinates  for  the  point 


I  **•<“( ■*£•<-!  +**4)/2 

Et,  )  yen  m(yei~i  +  ,Vf()/2  (54) 


j  yen-^yg'-.+yei 

2e»i  “  +  2t 


From  equations  (45)  and  (54),  it  is  possible  to  obtain  the  plane  of 
oscillation  *».•  froa  instant  l#_,  to  instant  t,  of  the  geometrical 
center  lines  if  the  j  *h  tub } .  Fuat  is, 


in  the  -“quations 


Aux  +  B„y  +  Cuz’*D,l 
A,,,B,,,Cu 


are  the  same  as  before. 


The  system  of  simultaneous  equations  ( 33)/ 44 ),(  55 ) 


,Atix+  Buy+Cti2=Du 
+  Buy  +  Cij*  =  Du 

^ A,ix  +  B,,y  +  C,iz**D,i 


makes  it  possible  to  solve  for  the  instantaneous  center  of  oscillation  of 
the  geometrical  center  of  a  jet  tube  .  This  is 

also  none  other  than  a  point  on  the  center  of  oscillation  of  the  jet 
tube.  It  is  possible  for  this  set  of  equations  to  be  solve!  through 

the  use  of  a  computer.  Therefore,  we  recognize  that  are 

already  known.  Later,  we  will  make  direct  use  of  them.  Fig.  8  is  an 
illustrative  diagram  of  the  1 nstantaneous  axis  of  rotation  of  the  jet 
tube  and  the  i nstantaneous  oscillation  center  of  its  geometrical  center 
line. 
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3,  Calculations  of  Moments  of  Force 

The  systems  of  moments  of  force  which  we  are  solving  for  here 
point,  under  differing  designated  wave  forms,  to  a  selected  series  of  ^ 
instants  (or,  one  might  say,  angles  of  oscillation)  to  figure  out  the 
overall  instantaneous  moment  of  force.  The  focus  of  this  article  is  on 
solving  for  the  moment  of  force-angle  of  oscillation  function 
relationship  (A/  =  /uJ))  .  Prom  this,  it  is  possible  to  obtain  the 

maximum  and  miniumum  (algebraic  values)  overall  moment  of  force  M ■•> 
and  Mmim  ,  as  well  as  its  aoment  of  force-angle  of  oscillation  function 
relationship.  Besides  this,  it  is  also  possible,  through  the  designation 
of  different  wave  forms  (normal  sine  waves,  sawtooth  waves,  and  square 
waves)  to  obtain  the  curve  of  the  functional  relationship  M .  As 
far  as  the  elaborate  separating  out  of  elastic  moments  of  force, 
asymmetrical  moments  of  force,  and  moments  of  friction,  and  so  on,  is 
» Moerned,  this  article  makes  no  additional  discussion. 

After  one  solves  for  the  instantaneous  center  of  oscillation,  the 
major  contradiction  in  a  solution  for  the  instantaneous  moment  of  force 
becomes  the  problems  of  the  point  of  action  and  the  direction  of  a  system 
of  forces  in  a  specified  space.  We  chose  to  make  use  of  the  two  methods 
of  Euler  transf ormations  and  four  dimensional  numerical  transforms,  and 
the  relative  difference  in  values  of  force  arms  which  we  obtain  in  the 
solutions  is  not  greater  than  0.8$.  Even  so,  we  only  introduce  one 
method  here. 

(l)  T'na  Print  of  Action  ijsed  to  Control  the  Direction  and  Amount  of 

Thrust 

From  Fig.  1  one  can  see  that  the  actuator  tube  movement  points  are 
~  ispaoV  vely  G,  and  G,  ,  the  points  of  action  for  motive  forces 
and  /ri  .  Moreover,  the  points  G,  and  G,  are  connected  to 
the  rigidity  provided  by  the  rod.  Therefore,  the  p  robin  b  n  >  a  n  one  of 
solving  for  the  spatial  coordinates  of  points  G,  and  G, 

(l)  Selection  of  a  dynamic  coordinate  system 
(a)  Translation  coordinate  system: 


VVT 


GX'Y'Z' 


- the  coordinate  axis  which  corresponds  to  the  fixed 

coordinate  system  QXYZ  is  always  maintained  parallel.  The  origin  point  G 
is  the  base  on  the  rod  center  line  for  point  G,  (or  . 


Fig.  9  The  Transf ormation  of  a  Dynamic  Coordinate  System  From  GX'Y 'Z'-GXWZ* 


(b)  Notary  Coordinate  System: 

GX'Y'Z '  -  the  angles  of  rotation  corresponding  to  the 

coordinate  system  GX'Y'Z '  are  and  their  rotational  order  is  as 

presented  below:  (See  Fig.  9) 

First,  we  rotate  through  angle  <p  around  axis  GX'  ,  obtaining 
coordinate  systme  GX,Y,  Z,  ■ 

We  also  rotate  through  angle  6  around  axis  GY ,  , 

obt a'  ning  the  coordinate  system  GX,Y,Z , 

finally,  we  rotate  through  angle  ip  around  axis  GZt  » 

obtaining  coordinate  system  GX'Y'Z'  .  The  angles  <p,< 5,y>  are 

called  Euler  angles. 


(2)  Euler  Ingles  and  Their  Transformation 


(a)  Solving  for  Euler  Angles 

In  the  calculations  or  Part  One,  we  already  obtained  the  result 


that  the  angles  of  rotation  around  the  various  axes  of  the  fixed 


coordinate  system  OXYZ  by  the  axis  line  of  the  jet  tube  are  p ,a,0. 

’Tow  we  mast  solve  for  the  Euler  ingles  which  correspond  to  each  instant. 

^ro'i  ^ig.  9  it  is  possible  to  obtain  the  transform  angular  velocity 
projection  relationship  for  the  GX'Y'Z'  coordinate  system. 


I  ♦ r. 


cosd'tinp 


cos &  cos W 


^ro «n  Table  1  we  get  the  following  series  of  relating  equations: 


dfi  d<p  ,  - 

"-'-■s*+s,nJ"3r 


da  d<5  .  .  dtf 

coy’  —  —r—  —  cos 93  •  - coeosini/i  -  -rf 

ar  at  at 


'  dd  dd  .  dti 

a,-  =  -7-  =  sin  gj  •  -7-  +  cos  6  •  cos  a  •  — ~- 

„  ai  at  at 


!d0  =  dgj  -f  sin  3  *  d  ^  ^ 

da  =  cos g>  •  dd— coed  •  sing?  ■  d^  ^  (56) 

d0*  singr>  •  dd  +  casd  •  cosgj  •  d$ 

door li ante  transformations  are  only  related  to  the  results  of 
rotational  <no  vements  and  do  not  consider  the  process  of  their  rotational 
aovements.  If  one  desired  to  solve  for  the  corresponding  angles  in  a 
certain  situation,  it  is  possible  to  recognize  that  one  always  begins  the 


rotational  movement  from  zero  degrees  arid  simul taneously  ou’ces  as »  o 


method  of  freezing  the  coefficients,  recognizing  the  use  of  appro ximatsly 
uniform  rotational  speed.  In  this  way,  equation  (56)  can  be  seen  as 
becoming  a  set  of  differential  equations  with  constant  coefficients  and  a 
zero  integration  constant,  for  purposes  of  integration.  Because  of  this, 


one  has: 


r0i  —  <T,  +  sin  d,  • 
a(  =  C06g?,  •  d,  —  coed,  sing?,  • 

0,  =  sing:,!?,  4-  coed,  •  cosg?,  •  V. 


(57) 


In  order  to  facilitate  the  carrying  out  of  iterative  substitution 
calculations,  we  change  these  to  the  form  shown  below: 


=  ai  •  sin<p,)/co6d, 

<P.  =  £.-*(  sin  <5, 

<3,  =  (<z,  +cos5<  ■  $in$\  •  i»,)/cos <p, 


(58) 


<p\  *=/J,id!  ==a,/coR95lt 

The  initial  value  for  iterative  substitution  can  be  taken  to  be: 

-  /\ 

•  finally,  we  solve  for  the  Euler  angles; 

(b)  Euler  Transf ormation 

The  purpose  of  the  transformation  is  to  carry  out  a  solution  for 
moments  in  a  fixed  coordinate  system.  On  the  basis  of  the  transformation 
order  which  is  given  in  Pig.  9,  we  obtain  the  results  shown  below: 

Transf  ormation  of  GX*Y*  Z'-*GX  ,Y,Z, 


—  sin  i> ,  0 

corv*  0  | 
/ 


(59) 


Transformation  of 


GA'.K.Z.-GX.V'.Z, 


coed  0  sin <5 

A,  =  [  0  l  0  j 

—  sin 5  0  coed 

Trans  for  nation  of  GX,Y,Z,^GX'Y'Z' 


( 60) 


“  0  cos  i p  —  sin  v 
L  0  sin  q>  cos  <p  - 


61  ) 


be 


ti>  tr  -msfo'-mation  natrix  from 
A , ,  Then ,  we  have 


GX'Y'Z ' 


to 


GA''K'Z' 


.  casiptoosd,  ,  — sintf,cosd,  , 

=  ,  sin )fifC06<p/  +  cos ^ / sin 3^  •  sin?),  ,  cos$,  •  oos<pt—  sin ^an<5,sin<)3t, 
sin^,  •  8in<p,  — oots sin d,  •  cos?,  ,  oosp(an^>i  +  sm?,  ■  sin<5,cos<p(> 


sin  <5, 

—  cos<3,  sin  <p( 

CO6<3,C06<p(  j 


The  transformation  of  point  G,,  in  coordinate  system  GX'Y'Z' 


r  xcu  *»  r  &G  i 

I  yen  I  =A,  0  | 

*cu  J  L  0  J 


The  transformation  of  point  Gu  in  coordinate  system  GX'Y'Z' 


■*ci»  "I  r  M 

y'cu  ~A,  Rc  ! 

.  J 


(j)  Tr  ins  fo’’  nat  ion  of  noting  rod  points  G,  and  G,  In  a  f  i  xed 
coordinate  system 

As  far  as  the  translation  operations  of  rigid  holies  are  ?  in  'em- 
the  a  in  ant  o  *  displacement  of  each  point  a  i  ?t  b  -*  -1  vni  vnl~n*  .  \  ■>  ■<  ,  •: 

>  »•  d  e  '*  t )  iol^e  for  the  amount  of  displacement  of  the  origin  pint  c 
.) "  the  translation  coordinate  system: 


Solution  of  the  s  haul taneoas  equations  (35)  and  (65)  allows  us  to 
obtain  coordinates  for  the  point  Q  : 

\  zg.-C-y  Vv ;  -  m  w 7,  jin ,  ^ 

S  XC(=>  K,2Ci+Q, 

t  _  T  _  ,  D 


in  the  equations 


vc,“  LiZc.+P, 

Ul  =  Kl+L]  +  \tVi  =  2lKi{Qt-x,i)  +  L/vPi-y.l)-z.J 

w,  =  (qi-xtiy  +  {Pt-ytly  +  zt\>wcl  =  ivl-H'c 


By  putting  equation  .63)  and  equation  (64)  together  respectively 
,-'.th  equation  (65),  it  is  possible  to  obtain  the  transfer  nation  of  point 
G,  (ju  in  tne  fixed  coordinate  sys tea. 


Xcu**x'cu+Xc, 
yCu’d’CM  +  Vc, 
ic„  =  ‘'h  +  'C, 


(67) 


*ti,  +  Xr, 

•>’«>• =  .vci(  +  y<.-, 
z  n  =  *c>,  +  m 


( 63) 


•t ion  of  the  forces  acting  to  control  the  amount  and 


'  -  i  on  of  forces  acting 

d  ; p i  i  'ate  the  axis  line  of  the  activator  tube. 
l  ■  .-di nates  for  the  dynamic  fulcrum  and  the  fixed 


ro"  their  directional  cosine. 

(69) 

4d”Xcu  Pot  Bh  —  y  c  i  itCii  =  rci,  —  H  o  ■  <  ^ 

•  =  x o it  * Bti  =  y Gtt  Po  =* zCn  —  H ou 1 


Arm  length  of  activator  tube: 


+  B, !  +  C,! 


Directional  cosine  of  forces  acting: 


i  70) 


,  cos  a ,  j  —  Ati  / Lmi 
F,,j  cos 0u  =  Bt,/Lxi 
oosyi  ( ™  Ct</ L„ 


f  co *a,,  =  A,JLy, 
also  F„<  COS0„  =  B,,/Ly, 

—  Cti  / Ly  , 


:  71 1 


(2)  Projections  of  the  forces  acting  on  the  various  coordinate 


axes : 


#F„,-F„C08alt 
|f,»i"F, jCos^i, 
Fi„  =  Fi(cosy,( 


also 


F,icasa,i 
F„  cos/?„ 
F„cosy1( 


7  ; 


(3)  Matrix  Calculations 

(l)  Solve  for  the  matrix  of  the  instantaneous  oscillation  :enter  o 
the  geometrical  center  line  of  the  jet  tube 

(a)  Radius  vector  of  the  i nstantaneous  oscillation  center  of  the 
forces  acting 

and  *  G,GU  or 


=  (jc6i(-xc,)/1000 

s=(>'cu-0'r,)/1000 

also 

^  rt*i  (.*<?!<  Jfri)/lO00 

•  j  r>y<~(ycn  —  yci)/ 1000 

=  (^Cl,  —  zc,)/l000 

r  (*G  ><  *C.  )/l000 

(b)  Overall  moment  of  force — Main  moment  and  i  i;s  relative  uimi 
'  ’  > j  ■;  j t  ions  of  the  main  moment  on  the  various  coordinate  axes  : 


\  *. 


w ..  *  y  r  r  ,F ,,,  ~r,,.Fr  -  l  r  ,,b  ,  4-  ,  r  f. 

t 

A/»,  ”  ^  (  r  ,,b  '  „f  ,,  )•(  r  „/•'  ,  y 


'  ..F  -r  .,/• 


•••  /■ 


•  v  W.  ♦  W  ♦  W. 


Vf,  -  Vf .  .  #, 


1  «  . 


A.,-B  ,C,.  -  r  -C  .  A.,  - 


ex*  [■_, 


_ A: f.M,. -  Af.,C_, 

-f  a/7:  i-  a/,:  %  +5.:  +  cv: 


V.'.v.-. 


7 ->  > 


jii*  j’M  '  l"1 )  and  (l~i  /  ,  we  can  then  solve  for  the  nooent  of 


f  tn*»  i  nstan taneoas  -ui  J  if  "it  it  lion 


M . ,  -  !  Mr ,  !ooip., 


•  i  !  < 


iff.  ip--  • .  f  .  a  j  aent  of  for  : 


■  <  .  i  i  o 


i  'l  e  in'  n  j  >  i  A  <  t  )  *■  h  )t  it  i  on 


j  •  .  ■  :  ’  it  i  m  -'enter  which  we  already  solved 

■1  I  i  ;  t  rlr'l  possible  tO  3  0 1 V  ‘ 


.  ■  ■  ;  I  •  I  I  '  I  l  i  i  <  !  >  e  "  )  t  1  t  i  )  0  : 


*u_  -  y-jn  - 


.  ■,  •  •  j  go  ’  n  rough  the  horizontal  pro  j-n't.  ion  din^r  in 
i  <  i  >'*  -  >Vt‘  :  >n  wni-h  appears  in  the  X —  Y  diagram, 

.  •  v  :  "  t  i  ■  .*  *nter  jf  oscillation  of  the 


o  c. ,  >i  s'  ,  v  *  so1  -tin  its  projection  equation: 


’  —  K  +  6, 


K_,-  B„,/A„,,b,  —  v,  K„,xc. 


•nter  line  >f  jet;  too  - 


♦  .  then: 


>*'  ".nn^in’  '  :  -  .'U  vi  i  mnye  .f  O-joO  degre  -•  s. 


■  '  i  ; 


B,  .  n  - "i ,  i  t  i  ?  ;j  o  ;ibu  :  to  precisely 


(3)  Calculation  of  Eccentric  Thrust 

Jet  tubes,  due  to  the  effects  of  additional  designated  oscillations 
ami  contained  pressures,  are  capable  of  causing  jet  tube  axis  lines  to 
give  rise  to  radial  displacement,  creating  thrust  eccentricity.  This  is 
an  important  technical  target  of  flexible  connector  heads. 

In  initial  configuration,  the  jet  tube  axis  line  equation  is; 

X-X._Y-Y._Z-Z.  ( -34) 

A.  B.  C. 


Tn  th  ?  ti nations: 

X,,Y,,Z,  - are  the  coordinates  of  a  certain  point  on  the  OZ 

a  x  l  s  . 

- is  the  directional  number  of  the  OZ  axis,  which  is 

always,  for  X,  =  Y,=*o,At— 5,=0,Z,  and  0,  equal  to  certain 

constants.  If  we  let  the  thrust  eccentricity  be  A  ,  and  we  solve 
e i  rat  ions  (35)  and  (34),  we  then  have: 


A’ 


X>t-X„  Y„  —  Yt,  Z,-Z, 

4,  Bt  C, 

A,  Bj  C, 


JA-B • 

V  A.B, 


B .  C, 

B,  C 


;ic.  a . 


C  A, 


_  \YttA,-X,lB,\ 
n/  B\  +  A] 


(35) 


:l;  »\ngala~  Velocity  of  Jet  Tube  Oscillation 

ly  Vc  f  ’  'i ->  '’unction  we  already  obtained  for  the  angles  of  rotation 


arouni  the  various  axes: 

0,-f(K  •  Ar)ia,-/(/C  •  bT)*6  =  f{K  •  A T)  K  =  0,\, n 


A T - equidist-.nt  stapling  intervals. 


cv.,  =  (0, —  AT 

«ax.  =  (a,-C(-i)/A7’  (36) 


[»] 


The  overall  angular  velocity  is: 


o,  -  Vo.*  +<», ;+<»,;  (87) 

(5)  Stroke  of  activator  tubes 

In  the  thirl  section,  we  already  solved  for  the  instantaneous 
lengths  of  activator  tubes  L.t  and  * L ,,  with  the  zero  position  length 
of  the  activator  tube  as  ~L,  .  Fro  a  this,  the  stroke  of  the  activator 
tubes  Is; 

r  ALm, *■  Lm,  —  Lt  (88) 

'  ALyi^ Lyi  —  L» 

(5)  Spatial  displacements  caused  by  the  effects  of  axial 
displacement  compensation  and  contained  jet  tube  pressures 

This  compensation  is  automatically  carried  out  when  the  command 
signal  is  zero.  After  that,  one  carries  out  measurements  and 
calculations  of  its  spatial  displacement.  This  is  also  another  important 
technical  parameter  of  flexible  connector  heads- 

With  jet  tubes  In  their  initial  configuration  (  pressure  zero, 
command  zero)  the  electrical  zero  value  of  the  mechanical  apparatus  and 
the  sensors  should  satisfy; 

After  carrying  out  compensation  for  the  effects  of  contained 
pressure,  this  should  satisfy: 

I  A"  t* —  y,i|<0.l  We  then  have: 

Axial  displacement: 

Radial  displacement: 


a z  =  z,;-z,:  } 

Ax=x,:-xt*t\  (89) 

a y-y.i-y,:  ] 


Note:  Units  la  all  the  equation  i  tni  salon  l.at  1  oa  s  ire  c  tho*n  biijw: 

Length:  The  r  ul  L 1  r,,.rn  a  i  »  a. it,  jn  u  •».  a  a ;  t, .  All  )Vi--  nils 
are  urn. 

Angles  and  angular  velocities:  these  respectively  use  degree;  and 
lig''vis  per  second  as  units. 

^  >  r  •  *  til  p.*  is  mi's:  these  respectively  use  kilogrvas  and  kilograms 
per  square  cm  as  units. 

Uonents  of  fore  -:  k  i  1  ogr  a  ns/  neter  . 

Time:  seconds. 
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